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FOREWORD 


Professor Chung-hua Wu pioneered the three-dimensional flow theory for turbomachines at Lewis 
Flight Propulsion Laboratory, NACA in 1950. He introduced the Sj and Sj families of relative stream 
surfaces and thus reduced three-dimensional flow problems to problems of iterating two solutions of two 
independent variables. The relaxation or direct matrix method was used for subsonic flows and the 
method of characteristics for supersonic flows. Laborious but accurate results were obtained without the 
benefit of modern digital computers. In the sixties Professor Wu developed a body-fitted, nonorthogonal 
curvilinear coordinate system to improve computational accuracy. In subsequent years Professor Wu and 
his colleagues at the Institute of Engineering Thermal Physics, the People’s Republic of China, developed 
shock-fitting and artificial compressibility methods for solutions in two- and three-dimensional transonic 
flows. Professor Wu’s theories were design tools used in aircraft engines such as the J69, JT-3D, Spey, 
RB211, JT9D, F404, etc. 

Since the early sixties Clemson University has been active in internal flow analysis. Through the 
support of the NASA Lewis Research Center in the early seventies, an inverse design method of the 
Griffith diffuser was developed. Initially the method was limited to potential flow. In subsequent years, 
the inverse method development at Clemson for internal flow has improved to include viscosity, 
compressibility and turbulence. Presently Clemson’s inverse solution method is used in design 
modification of the GE MS-7001F gas turbine using coal gas a s a fuel. 

Because of mutual interests in internal flows and ties to the NASA Lewis Research Center, I became 
familiar with Professor Wu’s work. In 1979, the year after U.S.A. and China resumed a normal 
relationship, I met Professor Wu in Beijing. He was the director at the Institute of Engineering Thermal 
Physics, The Academy of Sciences, People’s Republic of China, at that time. We discussed the 
possibilities of exchange visits and collaborations. 
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A project emerged with the following specific objectives: 

(1) To prepare a manuscript that summarizes the work of more than 100 journal articles on 

Sj and Sj methods and three-dimensional flow solutions in turbomachinery that Wu and his 
colleagues developed in the last 40 years 

(2) To give two lecture series on the above subjects, one at NASA Lewis Research Center and the 
other at the University of Cincinnati 

(3) To discuss, on a regular basis, research problems in internal flows with graduate students find 
faculty of the Department of Mechanical Engineering of Clemson University 

Thanks to the assistance of Dr. Melvin J. Hartmann, Director of University Programs at NASA 
Lewis, the above objectives materialized in 1990 with the support of a grant for NASA Lewis. Clemson 
University was privileged to have had Professor Wu and his wife, Professor Min-Hua Li, reside on 
Clemson’s campus from January 1990 until May 1990. The lecture series at NASA Lewis was held 
March 19-21, 1990, and April 16-17, 1990, in Cincinnati. The manuscript draft was completed prior to 
Professor Wu’s return to China. 

On the eve of publishing this report, I would like to recognize the efforts of Dr. Lonnie Reid, Chief 
of Internal Fluid Mechanics Division, and his colleagues who reviewed the manuscripts, and to express my 
appreciation to many others at NASA Lewis who made this report possible. 



It was a fulfilling experience and an inspiration for me to work with Professor Wu on this project. 
Persons who were involved with this project hope that this report will serve a useful purpose not only to 
document the work of computational fluid mechanics in turbomachinefy, but also to encourage those of 
us who continue to toil in turbomachinery research. 

On September 19, 1992, Professor Chung Hua Wu died in Beijing after a prolonged illness. While in 
the hospital, he read the final typing of the manuscript. We are sorry that he did not see his report 
released. The subject of this report is Professor Wu's lifetime effort. We hope this report will inspire 
those of us who toil in the field of turbomachinery. 


March 1991 


Project Coordinator 
Tah-teh Yang 

Professor of Mechanical Engineering 
Clemson University 
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covariant metric tensor of x 1 coordinate system 
contravariant metric tensor of x 1 coordinate system 
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enthalpy per unit mass of fluid, u + p/p 
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unit vector normal to stream surface 

cascade spacing (pitch) 

pressure 

any fluid quantity 

heat transfer to fluid per unit mass per unit time 
gas constant 
radius vector 

absolute cylindrical coordinates 
relative cylindrical coordinates 

relative stream surface passing through fluid particles lying on a circular arc upstream 
of or midway in blade row 

relative stream surface passing through fluid particles lying on a radial or curved line 

upstream of or midway in blade row 

entropy per unit mass of fluid 

absolute temperature 

time or circumferential thickness of blade 

blade velocity at radius r 

du = c y dT 

unit base vector and reciprocal vector 
absolute velocity of fluid 

angular momentum of fluid about axis of rotation 
relative velocity of fluid 

physical component of relative velocity tangent to x 1 
work done by fluid element per unit mass per unit time 
contravariant component of W 
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Wj covariant component of W 

x 1 general curvilinear coordinates (i = 1,2,3) 

z distance along turbomachine axis 

fi angle between W and its meridional component 

7 ratio of specific heats 

angle included by the coordinate lines x* and x^ 

P coefficient of viscosity, coefficient of artificial viscosity, or Mach angle 

S absolute vorticity, V x V 

r* viscous stress tensor 

p fluid density 

p artificial density 

a angle between z and I = tan"~ 1 (dr/dz) 

t normal, circumferential, or radial thickness of stream filament 

$ dissipation function 

9 stream function 

angular speed of blade 

d/ dx 1 partial differentiation of a flow variable on stream surface with respect to 
coordinate x 1 
Subscripts: 
c casing 

e exit station 

h hub 

i inlet 

I meridional component 

L.E. leading edge 


xi 



m 


mean (midchannel) 


n component in the direction normal to hub or casing 

p pressure surface of blade 

r,$?,z radial, circumferential, and axial component 

s suction surface of blade 

T.E. trailing edge 

Superscripts: 

° stagnation state 

* dimensionless quantity 


CHAPTER 1 


INTRODUCTION 

As a result of studying the effect of the radial equilibrium condition on the radial flow field in 
the axial-flow turbomachine (refs. 1 and 2), a general theory of three-dimensional flow in turbomachines 
(ref. 3) was proposed in 1950. It was intended for solving the three-dimensional flow in a turbomachine 

(1) Having arbitrary hub and casing shapes — the theory applicable to axial-flow, radial-flow, 
and mixed-flow turbomachines 

(2) With a finite number of blades which have finite thickness and arbitrary shape 

(3) With fluid moving through it at a high speed — the speed of flow being purely subsonic or 
supersonic 

The fluid flow through the stator and rotor blade row was assumed to be steady with respect to the 
stationary blades and rotating blades, respectively. It was proposed to obtain steady flow relative to the 
blades by an iterative solution between two families of relative stream surfaces. The families were the 
family and the S 2 family. The problem of determining the flow field with three independent variables 
was reduced to a number of flow fields having only two independent variables. Thus, the purely subsonic 
or purely supersonic flow along an S t or S 2 relative stream surface could be accurately solved by the 

techniques available at that time. The relaxation or direct matrix method was used for the 
subsonic flow, and the method of characteristics was used for supersonic flows. References 4 to 10 
contain the solutions obtained by these methods and the approximate solutions obtained by series 
expansion in the circumferential direction from a mean streamline on the Sj stream surface in axial flow 
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turbomachines and centrifugal compressors. Almost all of these calculations were obtained by 
mechanical, digital computers. This was a laborious endeavor and consumed a great deal of time, but 
provided insight into the characteristics of three-dimensional compressible flow in turbomachines. Those 
insights were useful to the development of the three-dimensional flow in turbomachine theory. 

The flow computations reported in references 10 and 11 are probably the First two turbomachine 
calculations ever performed on large-scale, high-speed, electronic digital computers. Reference 10 includes 
calculations along the mid-channel S 2 relative stream surface for high subsonic flow. Reference 8 
contains computations along the mid-span Sj relative stream surface for incompressible flow at design 
and off-design inlet angles. The number of interior grid points for the two problems was 400 and 200, 
respectively. The fourth degree differential formula was used. It took approximately 60 hr on an IBM 
CPEC computer to factorize the coefficient matrix into lower and upper triangular matrices for the 
turbine problem. The compressor problem was also factorized in 60 hr on an IBM 604 computer. The 
gas turbine problem was also solved on an UNIVAC computer later and took a relatively shorter period 
of time — 11 min for factorization and 2.5 min for each cycle of stream function calculation. 

With the advent of the faster, modern digital computers, solutions for the subsonic flow along the 
and S 2 stream surfaces were obtained by many turbomachine investigators (refs. 12 to 26). 

Solutions obtained on IBM 360, 370, KDF9, and Facom 230-26 took about 0.5 to 35 min. In addition to 
the S x and S 2 stream surface flow solutions, quasi- (refs. 22 to 25) and full-three-dimensional flows 
(ref. 25) were obtained on IBM, Facon, and CDC 7600 in about 2 to 16 min (refs. 22 to 25). The 
surfaces were assumed to be surfaces of revolution in calculating the quasi-three-dimensional flow. 

Along with the development of high speed digital computers, the development of mathematical 
calculation techniques continued. One major development of the latter was avoiding the inconvenience 
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and inaccuracy caused by unequal grid spacing near an arbitrarily shaped boundary surface, which exists 
in all practical applications. One approach used more accurate, high-order differentiation at such points 
(ref. 27). A second approach used body-fitted finite elements to discretize the differential equations (refs, 
19, 20 and 25). A third approach used the body-fitted, nonorthogonal, curvilinear coordinate system 
(refs. 28 to 33). A second major development in mathematical techniques was the solution of transonic 
turbomachine flow using the time-marching method (refs. 34 to 36), separate-region .calculation and shock 
fitting (refs. 37 to 39), and artificial compressibility (refs. 40 to 42). 

This report is a slightly expanded version of a series of lectures given at NASA Lewis and at 
Cincinnati University during the spring of 1990. 

Chapter 2 briefly reviews the fundamental aerothermodynamic equations governing three- 
dimensional flow in turbomachines. The equations, in the beginning, are for the most general case of 
unsteady flow of a viscous fluid relative to a rotating blade. The two independent thermodynamic 
properties selected were entropy and relative stagnation rothalpy of the fluid. The latter is a 
thermodynamic property introduced especially for calculating three-dimensional flow in the rotating blade 
row. The implication of assuming steady absolute flow in the stator and the steady relative flow in the 
rotor at the same time is also discussed. After a discussion on the effect of viscosity in the governing 
equations, a practical method for considering the viscous effect on the flow field is given. 

Chapter 3 presents the basic idea of expressing flow variables on a general stream surface in terms of 
two independent variables, i.e., the two coordinates. This chapter continues conservation of mass in a 
fluid element. In such a case, the governing equations naturally contain a stream filament thickness 
term, a general function of the two coordinates. The flow along the stream filament is obtained by 
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solving the principal equation with the stream function as the single dependent variable. Procedures for 
solving direct (analysis) and inverse (design) problems are described. 

To have high accuracy in the finite-difference approximation at grid points near a curved boundary 
wall and to satisfy the boundary condition at the curved wall accurately, the body Fitted, general, 
nonorthogonal, curvilinear coordinate system is used for both Sj and S 2 stream filaments. By using 
tensor calculus — the continuity equation, the vorticity equation, the dynamic equation, the energy 
equation, and the principal equation are easily expressed in terms of the general, nonorthogonal, 
curvilinear coordinates and corresponding nonorthogonal velocity components. These equations and 
methods of solution for flow along Sj and S 2 stream filaments are given in Chapter 4. 

The first part of Chapter 5 presents a simple, approximate solution for subsonic flow along the Sj 
stream filament of revolution by circumferentially extending the known solution on the mid-channel 
streamline. Also in this chapter a simple, approximate three-dimensional solution is obtained by 
circumferentially extending the known values on the mid-channel stream surface. The second part of this 
chapter presents results obtained in subsonic Sj solutions employing H type and C type 
nonorthogonal, curvilinear coordinates. The third part describes procedures for quasi-three-dimensional 
blade design and full-three-dimensional analysis of given blades. A comparison of the calculated three- 
dimensional flow field and measured data is also included. 

Chapter 6 describes several relatively quick methods for calculating the transonic flow along S t and 
S 2 relative stream filaments. The method of separate-region calculation with shock-fitting, elliptic 
solution of the stream-function principal equation, to which artificial viscosity is introduced in the density 
term to stabilize the transonic calculation for both the Sj and S 2 stream filament, and the elliptic 
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algorithm for the inverse solution of S 2 flow (V^r prescribed), which is modified for obtaining a sharp 
shock discontinuity, are presented. The calculated results are compared with experimental data. 

Applying the quick solution methods described in Chapter 6, the quasi-three- and full-three- 
dimensional transonic flow solutions in two compressor rotors were obtained rad are presented in 
Chapter 7. The solutions are presented with emphasis on the convergence process and the geometry of 
individual S 1 and S 2 stream filaments obtained in the three-dimensional solution. These solutions are 
also compared to experimental data and are included in Chapter 7. 

Based on the analytical solutions of three-dimensional subsonic and transonic flows and their 
respective experimental data, practical methods for three-dimensional turbomachine blade design and 
blade element test data correction emerged and are proposed in Chapters 6 and 7. 
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CHAPTER 2 

FUNDAMENTAL AEROTHERMODYNAMIC EQUATIONS GOVERNING THE 
THREE-DIMENSIONAL FLOW IN TURBOMACHINES 

2.1 Basic Aero thermodynamic Equations Governing the Three-Dimensional 
Flow of a Viscous Fluid Through a Stationary Blade Row 

The general basic aerothermodynamic equations governing the flow of a viscous fluid through a 
stationary blade row which were formulated in reference 1 in connection with the calculation of a radial- 
equilibrium condition for the design of turbomachine blades are as follows: 

Continuity Equation: From the principle of conservation of matter, the equation of continuity is 

^ + V(/»V)= 0 (2.1) 

<Jfc 


or 


r.v + £L' = o 

Dt 


(2.1a) 


Dynamic Equation: Newton’s second law of motion is expressed for viscous fluid by the Navier-Stoke’s 
equation. 



-Vp + n 


V 2 V + 2 V (V • V) + 
3 


( 2 . 2 ) 


where . . . represents high order terms due to viscosity change with temperature. 
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Energy Equation: From the first law of thermodynamics, one form of the energy equation is 


+ p D ( 1 /P) 

Dt Dt 



in which 


du = c v dT 


(2.3) 


(2.4) 


q = I V ■ (kVT) 

P 


(2.5) 


and 


9 = p 


2V • (V- V)V + (Vx V) 2 - 2(V ■ V)(V • V) - * (V • V) 2 


( 2 . 6 ) 


For turbomachine calculations it is convenient to take specific stagnation enthalpy H and specific 
entropy s as the two independent thermodynamic properties defining the thermodynamic state of the 
gas. They are related to other thermodynamic properties and velocity by 


H = h + lv 2 = (u + p/p) + I V 2 
2 2 


(2.7) 


and 

T ds = du + p d(l /p) (2.8) 

By the use of equations (2.2), (2.7), and (2.8) the following forms of dynamic equation and energy 
equation were obtained (ref. 1): 
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(2.3a) 


DH 

= q+I^ + i + V- 

P 

y 

v 2 v + iv(v • v) 

+ • • • 

Dt 

p dt p 

P 

3 J 

, 


In the case of steady invicid flow, equation (2.2a) becomes 

V x (V x V) = VH - TVs ( 2 - 2b ) 

which is the Coroco equation originally deduced for the investigation of flow with shock and vorticity. 

The second law of thermodynamics states 


T > q (2 9) 

Dt 


Combining equation (2.8) with the first law of thermodynamics, equation (2.3), yields 


T 


Ds 

Dt 



which conforms with the second law of thermodynamics, equation (2.9). 


( 2 . 10 ) 
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From the preceding general, basic equations the following can be noticed: 

(1) In the stationary frame of reference the flow unsteadiness is represented by the partial derivative, 
with respect to time, of density in the continuity equation, of velocity in the dynamic equation, and of 
pressure in the energy equation. 

(2) Stagnation enthalpy is, in general, affected by the viscosity of the fluid through the last two 
terms on the right side of equation (2.3a). If the last term on the right side of equation (2.1a), the 
viscous force per unit mass of fluid, is denoted by Fj, then the two viscous terms in equation (2.3a) are 
*/P and V-F r Thus the effect of viscosity on the stagnation enthalpy is not represented by V F f 
alone. 

2.2 Effect of Viscosity on Basic Equations 

In order to more clearly see the effects of viscosity on the changes in stagnation enthalpy and 
entropy, the dynamic equation (2.2a), the energy equation (2.3a), and the entropy equation (2.10) are 
further examined using the stress tensor r (ref. 43). 

Newton’s Second Law of Motion 

Let the resultant stress (force per unit area) acting on the surfaces of an infinitesimally small fluid 
element be denoted by x (fig. 2.1) where r is related to the hydraulic pressure p, uniform in all 
directions, and a viscous stress r' by the following equation: 
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( 2 . 11 ) 


or 




r xy 

h" 1 

hi 


°k 

r xy 

r xa 


-P 

0 

0 

1C = 

r y* 

°y 

r y* 


r yx 

°y 

r y* 

+ 

0 

-p 

0 


r 

T „ 




T 

O' 


0 

0 

-p 


' ax 

*y 

a 


■X 

*y 

2 






*n = - pS a + r ij 


( 2 . 12 ) 


where 


s V* i=j 

IJ 0, i 1* j 


The vector and tensor form of Newton’s second law of motion are, respectively 


DV 

Dt 


-iVp + 1 V V 

P P 


(2.13a) 


and 


Dt dt j dxj p dxj p dxj 


(2.13b) 


where Einstein’s summation convention is used. 

By using equations (2.1), (2.2), and (2.3), Newton’s second law of motion can also be put into the 
following form: 
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(2.14) 


— - V x (V x V) = -VH + T x Vs + 2vV 
fit p 

First Law of Thermodynamics 

The first law of thermodynamics is, in general, expressed by 



w 


(2.15) 


The time rate of work done by the force acting on the fluid element surface, per unit mass of fluid, 
as seen by a stationary observer is 


* _ » a M _ _ i a H ^ » a h v il 

P dxj p dxj p dxj 


or 


where 


w = Iv.(jr-v) = -2v(pV) + 2 V ■ (ir'.v) 

P P P 


x'.V = r&V, 


The first term on the right side of equation (2.16b) can be written as 


(2.16a) 


(2.16b) 


(2.17) 


2 V-(pV) = 2 (V-Vp + pV V) = 2 

P P P 


Dp 

,Dt 


dp 

dt) 


+ E v-v 

p 


li 



Substituting the continuity equation (2.1a) into the preceding equation results in 


-V.(pV) 

p 


if Dp _ dp _ p 1 D p _ D(p / p) _ 1 dp 
p Dt dt j p p Dt Dt p d t 


Substituting the preceding equation into equation (2.16b) results in 


w = D ( p / P A - + I V • (x'-V) (2.18) 

Dt p dt p 


To this observer the time rate of increase of internal energy is 


De _ D u 

Dt Dt [ 2 ; 


(2.19) 


Substituting equations (2.18) and (2.19) into equation (2.15) yields 


22 = 12 ? + q + IV • (x'V) . (2.20) 

Dt p dt p 


The energy equation (2.20) can be put into a slightly different form by expanding the last term on 
the right side of the equation a s follows: 


i 

_ 1 

dV: 

*' i{ l + v 5 

dir'y 

p dxj 

P 

** ax, ‘ 

dx£ 


( 2 . 21 ) 
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But 


, dV . 

’“s = t 


i _ 
P 9xj 


Therefore, the energy equation (2.20) can be written as 


= + q + i + Ff -v 

Dt p dt p 


which is the same as equation (2.39). 

The energy equation (2.20) can also be transformed to the form of equation (2.3) as follows by 
multiplying terms on the left side of equation (2.13b) by Vj: 


Dt 2 1 dxj 1 dxj 


Using equations (2.21), (2.22a), and the continuity equation (2.1a), equation (2.14) is transformed into 
the following equation: 


* = p D (Vg) _ ♦ _ J5 V? 

[ Dt p\ Dt 2 
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Substituting equations (2.25) and (2.19) into the first law of thermodynamics, equation (2.14) results in 


Du , D(l/p) . , * 

— + p \ ill - q + _ 

Dt Dt p 

which is identical to equation (2.3). 

The physical meaning of this form of the energy equation can best be seen by deriving the equation from 
the point of view of an observer moving with the fluid element. To this observer the time rate of the 
increase of its internal energy is 


De _ Du 
Dt Dt 


(2.26) 


and the rate of doing work against the surroundings by 


unit mass of the fluid element per unit time is 


w 


( 


1 3Vi 

/ ij 




_ r P(V/>) _ £ 

Dt p 


(2.27) 


In equation (27) is the rate of work done by unit mass of fluid element as seen by an observer 
moving with the fluid element and are the coordinates moving with the fluid element. 



Substituting equations (2.26) and (2.27) into equation (2.15) yields 

(2.28) 

Equation (2.28) is exactly the same as equation (2.21) which is obtained by modifying the energy 
equation (2.26) obtained from the point of view of a stationary observer with the use of Newton’s second 
law of motion (eq. (2.13)). 

It is important to notice the following: 

(1) To a stationary observer, the rate of work done by unit mass of fluid against viscous forces and 
the rate of increase of internal energy are, respectively, 


Du . I D(1 /,) _ * 


Dt q 


• V) and D (u + V V 2 ) (2.29) 

P Dt 

These two terms appear in the corresponding energy equation (2.20). 

(2) To an observer moving with the fluid, the rate of work done by the unit mass of fluid against 
viscous forces and the rate of increase of internal energy are, respectively, 



$ 

P 


and 


Du 

Dt 


These two terms appear in the corresponding energy equation (2.28). 
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In the derivation of the energy equation from the first law of thermodynamics, the expression for 
work done by the fluid element against its surroundings and the increase in the internal energy of the 
fluid element must be written for the same observer, either stationary or moving with the fluid. If one is 
written for the stationary observer and the other is written for the observer moving with the fluid, the 
resulting equation is erroneous. Unfortunately this kind of mixup has appeared in some publications. 

Two-Dimensional Laminar Boundary Layer Flow In order to clearly see the effects of fluid 
viscosity and heat transfer on the energy equation and entropy equation, the following general analysis of 
steady laminar boundary layer flow is made. The velocity and temperature distributions in the boundary 
layer and an infinitesimally small element whose width and height are Sx and 6y are shown in 
figure 2.2. The x and y coordinates are chosen, respectively, to lie along and perpendicular to the 
tangent of the blade surface. The approximate relation commonly used to treat boundary layer problems 
is employed in the analysis. Fluid pressure may have a gradient in the x direction. 

(1) To a stationary observer, the time rate of work done by the fluid, per unit mass, is 


' v ,) + + pV y) 

oy 

v-U'-v)] 

and after employing the continuity equation 

w - p (p/*> - Iv .(y.v) 

Dt p 



(2.30) 


(2.30a) 
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The time rate of heat transfer to the fluid element, per unit mass, is 


q = 1 V 
P dy 


k£E 

l dy) 


(2.31) 


Substituting equations (2.30a) and (2.31) into the first law of thermodynamics equation (2.15), the 
following is obtained: 


1 d 

k £l| - 

D(p/ p) _ 1 

3 K) 4 V «) 


D 

f 

. v a 
11 •+■ 

P dy 

dy J 

Dt p 

dx dy 


Dt 

u 1 — — 

2 J 


Rearranging terms gives 


DH 

Dt 


d 

r \ 

k dh 

+ 

a K) + a K) 

dy 

.Cp dy, 


d* dy 


(2.32) 


(2.33) 


In equation (2.33), 
then becomes 


3(rV,) 

dx 


is much smaller than 


d(rV x ) 

and can be neglected. Equation (2.33) 

dy 


dH 

~dt 


1 

d\ k ah 

, d 

f w 

av_ 

— 


+ _ 

**v x — - 

P 

dylcp dy ; 

dy 

dy JJ 


(2.34) 
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If n = constant, Pr = 1, and the difference between V and V x is neglected, then equation (2.34) 
becomes 


dH 

It 


ft — + — — 


!t — 
p ay 2 


(2.35) 


H = constant everywhere in the boundary layer is a particularly useful solution to the above equation. In 
this case, the wall temperature T w equals the fluid stagnation temperature T 0 . When c p of the gas is 


a constant. 


T = T + = constant 

2c P 


(2.36) 


Differentiating equation (2.36) results in 


dT _ v av 

dy c p dy 


It can be seen from equations (2.37) and (2.31) that 


(2.37) 


(a) at the lower boundary of the boundary layer, V = 0 and dV/dy > 0, therefore, r > 0 and 

arr/dy = q = o ; 

(b) at the upper boundary of the boundary layer, V = Vj and dV/dy — 0, therefore, r = 0 and 

ar/ay = q = o ; 

(c) between these two boundaries, dW/dy > 0, r > 0, dT/dy < 0, and therefore, heat transfer is in 
the y direction. 
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Since H = constant, equation (2.34) indicates that, under the present approximation, the amount of 
heat transferred into (or out of) the fluid element is equal to the work done by the fluid element on the 
surrounding (or by the external viscous shear stress on the fluid element), i.e., 


B 

r i 

k^I 

= - -K) = - 

f \ 

d^ 

T 

-k- 

dy 

, dy j 

dy' 

dy J 



(2.38) 


Of the two terms on the right side of equation (2.38) the first one is always negative and the second one 
is always positive. 


2. To an observer moving with the fluid the time rate of work done by the fluid element, per unit 
mass, on its surrounding is 


w = pEiiM-: 

Dt p 

where £ and rj are coordinates moving with the fluid, with the velocity W. Since 

dW v /d{«dW f/dt) 




dW. 


— i + — '1 
dr, B£ 


~ r D(Vrt _ r aV * 

Dt o Bv 


(2.39) 


$ 



dt , 




> 0 


(2.40) 


Substituting equations (2.39) and (2.32) into the first law of thermodynamics, the following is obtained: 
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(2.41) 


Du 

Dt 


1 d 
P dy 


( k - 


p D(l/,) _ r W,' 

l dy) 


„ Dt p dy j 


Substituting equation (2.41) into the entropy equation (2.61) results in 


ip Ds 

1 

d 

, 3T 

X 




k — 

+ r 

Dt 

P 

dy 

. dy, 

dy J 


(2.42) 


The preceding two equations give, respectively, rate of increase of internal energy and entropy in the fluid 
element due to the work done on the fluid element by viscous force acting on the fluid element. 

In the particular solution above, substituting equation (2.38) into equations (2.41) and (2.42) results 
in the following equations, respectively: 


2 

Du = l v ^_p EKVf) = ft - p (2-43) 

Dt p X dy Dt p dy 2 2 Dt 


and 


v 2 

T Ds = _ 1 v = (2.44) 

Dt p X dy P dy 2 2 

The first term on the right side of these two equations is always positive because dr/dy is always 
negative. 

In this example it is quite clear that 
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(a) the rate of work done by the viscous forces on the fluid element, per unit mass, seen by a 
stationary observer at rest is 1/p d(rV x )/dy = 1/p {rdV x /dy + V x dr/dy), where the first term is 


always positive and the second term is always negative. 

In the particular solution above, because the work done by the fluid element and heat transfer to the 
fluid element just cancel each other, the stagnation enthalpy remains constant along the relative 
streamline. 

(b) The rate of work done on unit mass of fluid by the viscous force seen by an observer moving 

with the gas is r(l/p) dV x /dy, which is always greater than zero. 

In the particular solution above, (q — w) = — 1 V v — — in which the first term is 

P dy Dt 

always greater than zero, thereby causing Du/Dt and Ds/Dt always to be greater than zero. 

It is quite obvious that the work done by the fluid element seen by a stationary observer and by an 
observer moving with the fluid element are not the same. In writing out the energy equation, they should 
be used with the increase in internal energy of the fluid element seen by the same observer. 




2.3 Basic Aerothermodynamic Equations Governing the Three- Dimensioned Flow of a 
Viscous Fluid Through a Rotating Blade Row 

General basic aerothermodynamic equations governing the fluid flow through a blade row rotating at 
a constant angular velocity were formulated in reference 3 for a nonviscous fluid, in which the entropy s 
and a new thermodynamic property I, first called “modified total enthalpy” and later named “relative 
stagnation rothalpy” (ref. 44), were taken as the two independent thermodynamic properties defining the 
thermodynamic state of the gas. Later this formulation was extended to viscous gases in reference 43. 


Continuity Equation: 


From the principle of conservation of matter, the equation of continuity is 
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^ ! + v-(pW) = o 

dt 


(2.45) 


or 


V • 


W + 


D In p 
Dt 


= 0 


Dynamic Equation: For a blade rotating at a constant angular velocity u about the z 
second law of motion is 


5?! = w J r + 2 «xW = -ivp + iv-f' 
Dt p p 


Where *■' is the viscous stress tensor acting on a fluid element. 


The relative acceleration in equation (2.46) can be written as 


DW = £?[ + (W-V)W = aW + IvW 2 - W x (V x W) 
Dt dt dt 2 


Rothalpy and relative stagnation rothalpy are defined, respectively, by (ref. 44) 


i = h - 


M 2 

2 
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(2.45a) 
is, Newton’s 


(2.46) 


(2.47) 


(2.48) 



(2.49) 



(wr) 2 W 2 

2 T 


By using equations (2.4), (2.8), and (2.49) the following form of Newton’s second law of motion is 
obtained: 


— - W x (V x V) = - VI + TVs - Ivt' (2.50) 

dt p 

Energy Equation: The energy equation for fluid flow passing through a rotating blade can be obtained 
from the first law of thermodynamics in the same manner as in the case of fluid flow passing through a 
stationary blade row. First, from a stationary observer’s point of view, the rate of work done by the fluid 
element, per unit mass, against its surrounding is 

w = - EF • V = - EF • (W + U) 

= - iv • (x*W) — EFawr 

P * ( 2 . 51 ) 

= Iv • (pW) - Iv • (i'W) - w D ( V,r ) 
p p Dt 


where U is the blade speed wr. To the same observer the rate of increase of the fluid element’s internal 
energy per unit mass is 


De _ 

D 

f \ 

a. V 2 

u -|~ 

D 

W 2 + 2 Wawr + <A 2 
11 + W 

Dt 

Dt 

2 J 

Dt 

2 


(2.52) 
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Substituting equations (2.51) and (2.52) into equation (2.15) yields 

= + q+lv-(»'W) . (2.53) 

ut p at p 

Just as in the case of fluid flow passing through a stationary blade row, the energy equation in the 
form of the preceding equation can be transformed into the form similar to equation (2.28). For instance 
the dot product of W and equation (2.46) is 


D W 2 

Dt ~2 


W • (<A - 2« x w) = W, - W; 

dxj dxj 


(2.54) 


The last term in the preceding equation can be written as 


W, 


■*« _ 4'« w i) _ aw, 


dx t 


8xj 


*1 


(2.55) 


and the dissipation function for the fluid element moving with respect to the rotating blade is 


3Wj 


(2.56) 
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By using equations (2.55), (2.56), and the continuity equation (2.45); equation (2.51) can be transformed 
into the following equation: 


w = 

p D(i/<o _ r 

D W 2 -w 2 r 2 _ m 


i Dt p , 

Dt 2 Dt 


(2.57) 


Substituting equation (2.57) and (2.52) into equation (2.14) yields the following form of the energy 
equation: 



r du/p) 

. Dt p ; 


(2.58) 


Again the physical meaning of this form of energy equation can be seen more clearly by deriving the 
equation from the point of view of an observer moving with the fluid element. To this observer the time 
rate of the increase in the internal energy is given by equation (2.26) and the time rate of work done by 
the fluid element against its surroundings is 


1 aW i D(1 Ip) *' 

W = — _ Wtz = n v ' — 

P ' 


Dt 


(2.59) 


where £ is the coordinate moving with the fluid. 


Substituting equations (2.59) and (2.26) into equation (2.14) yields equation (2.58). 
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Entropy Equation: From the second law of thermodynamics, 


T 


Ds 

Dt 


> q 


(2.60) 


substituting (2.58) into the entropy equation (2.8) the following is obtained: 



(2.61) 


Equations (2.53) and (2.61) are two important equations. They appear quite similar to equation 
(2.20) and (2.10), respectively, but the total derivative with respect to time now means following motion 
along the relative streamline, and the partial derivative with respect to time now refers to the derivative 
at a coordinate point relative to the rotating blade. 


2.4 Steady Fluid Flow Through a Stationary Blade Row and 
Steady Relative Flow Through a Rotating Blade Row 

Under steady operating conditions, fluid flow through a single stationary blade row is steady and the 
unsteady terms in the governing equations can be neglected. Similarly under steady operating conditions 
fluid flow through a single rotating blade row is steady and the unsteady terms in the governing 
equations can be neglected. However, even in a single-stage turbomachine, there is always a stationary 
blade row upstream (in the case of a turbine) or downstream (in the case of a compressor) of a rotating 
blade row. Usually they are spaced not too far apart and the fluid flow relative to either blade row is 
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unsteady. Because of the mathematical difficulty, in practically all of the design calculations and analysis 
calculations, fluid flow relative to either the stator or rotor, is assumed to be steady. 

Steady absolute flow in stationary blade rows means that the unsteady terms in governing equations 
(2.1), (2.14), (2.20), and (2.28) at a fixed coordinate point, say (r,0,z) Q are equal to zero. 

In the case of fluid flow through a rotating blade row, steady relative flow means that all of the 
partial derivatives with respect to time in governing equations (2.45), (2.50), (2.53) and (2.61) at a fixed 
coordinate point, say (r,^,z) 0 , with the origin of the coordinates fixed in the blade (i.e., <f> = $ -wt), are 
equal to zero. Furthermore, since the absolute velocity V is related to the relative velocity W by 

V = W + U (2.62) 


or 


V r = W r 

v, = w, 
v, = W, + «r 


(2.63) 


When the relative velocity at a relative coordinate point (r,^,z) is steady, the absolute velocity at a 
relative point, V r ^ B is also steady with respect to the rotating blade, so is V x V. 


Because the absolute flow is calculated for the stationary blade row, whereas the relative flow is 
calculated for the rotating blade row, there is an abrupt change in the tangential component of the fluid 
velocity, and consequently, in the streamline when the fluid motion is referred to the different coordinate 
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system moving from one blade row to the next (fig- 2-3). However, the projection of the streamlines on 
the meridional plane is continuous because the meridional components of the absolute and relative 
velocities are continuous (fig. 2.4). 

In the following presentation, as well as in computer codes, only the governing equations for steady 
relative flow through the rotating blade row will be given. It is understood that when the blade row is 
stationary, u — 0, <p — ► 0, W — ♦ V, I — ► H, and $'—►$. 

2.5 Viscous Terms in the Governing Equations 

An analysis will now be made on the magnitude of the viscous terms in the governing equations. 
Continuity Equation: For steady flow the continuity equation (2.45) becomes 


V- (*>W) = 0 ( 2 - 64 ) 

Equation (2.64) does not contain a viscous term. The effect of viscosity on the fluid flow comes through 
the entropy increase in the flow in the following equation of fluid density: 


— 

P* 


(i — W 2 /2 + U 2 /2)b 


(i - W 2 /2 + U 2 /2)„ 


7-1 


s b -s. 

. R 


(2.65) 


The effect of entropy increase in density is quite large and consequently cannot be ignored. 

Energy Equation: Of the two forms of energy equation (2.53) and (2.58), it is more convenient to use 
(2.53). For steady flow it becomes 
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( 2 . 66 ) 


El * q + 1 V • (»'• W) 
Dt p 


In the core region the viscous stress and heat transfer are negligible. In the boundary layer region near 
the blade surfaces and hub and casing walls, if the boundary layer is laminar, the Pr number of the fluid 
is equal to unity, and the boundary walls are adiabatic walls, the viscous work term and heat transfer 
term cancel each other. In actual turbomachines the boundary layer flow is turbulent and the Pr number 
is different from 1, and the summation of these two terms will not be equal to zero, but the magnitude is 
expected to be small. The following equation is usually considered to be a good approximation for the 
entire flow region: 


DI 

Dt 


= 0 


Dynamic Equation: For steady flow equation (2.34) becomes 


(2.67) 


W x (V x V) = VI - TVs + - Vx' (2.68) 

P 

When I is taken to be constant on all streamlines, the magnitude of VI depends on the magnitude 
of VI at the inlet. Vs is quite small in the core region, but quite large near the solid wall. For 
instance the radial entropy profile at the exit of a rotor blade row may look something like that shown in 
figure 2.5. From the dynamic equation in the radial direction. 
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When the viscous stresses are neglected, equation (2,68) is simplified to 


W x (V x V) = VI - TVs 
Equation (2.69) in the radial direction is 


W 


% 


dW 




dr 




fW + Wi f^ + f»_ T * 

r dr di dr dr 


(2.69) 


(2.70) 


The effect of radial entropy gradient on the radial variation of velocity is indicated in figure 2.5. 

The inclusion of entropy gradient in the dynamic equations is necessary for a better prediction of velocity 
variations near the solid wall. 

Entropy Increase: Evaluation of entropic increase along the streamlines by equation (2.61) requires 
a solution of the complete set of governing equations for viscous fluid. At the present time approximate 
value of entropy increase along the streamlines may be estimated by an appropriate value of the 
poly tropic exponent of a poly tropic process which represents the actual flow process (ref. 1). It may also 
be estimated by the pressure recovery factor and the isentropic rotor efficiency obtained in experimental 
investigation and tests (ref. 47). For the flow through a stator blade row the entropy increase across the 
blade row is calculated from the recovery factor in stagnation pressure as follows: 



(2.71) 


or 


30 


1 II 



(2.71a) 


s — s b = R 2_ In . — 

7-1 P c . 

For the flow passing a rotor blade row, the increase of entropy across the blade row is calculated from the 
isotropic rotor efficiency as follows: 


8 b ~ 8 a = R 


. In 


T b -/T. 


7-1 1 + f(T b ./T,. - l) 


(2.72) 


where 


£ = rf s for compressor rotor 
=l/i/| 8 for turbine rotor 


2.6 Some Remarks on the Energy Equation and the Entropy Equation 

The First Law of Thermodynamics and the Second Law of Thermodynamics are two important 
physical laws governing the flow of a compressible viscous fluid in a turbomachine. The following 
remarks are made here regarding the energy equation and entropy equation derived from these two laws. 

1. Energy equations in two different forms were described in Section 2.1 and 2.2 for flow through, a 
stationary blade row and a rotating blade row respectively. In one form of the energy equation, equation 
(2.28) or (2.58), the increase of internal energy of a fluid element along the streamline is given directly by 
the heat transfer from the surroundings to the fluid element minus the work done by the fluid element 
against its surroundings. The equation is the same as the one universally used in thermodynamic 
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calculations, except that there is an additional work done term by the fluid element against the external 
viscous force. In the other form of energy equation, equation (2.2) or (2.53), however, the increase of 
stagnation enthalpy or stagnation rothalpy is given by three terms, namely, the heat transfer to the fluid 
element, the work done by the fluid element against the viscous forces, and an insteady pressure term. It 
shold be noted here that the work done by the fluid element against external pressure is included in the 
first form of the energy equation, but not in the second form of the energy equation. (Equation (2.18) 
will help to explain this difference.) 

2. It should be emphasized that the heat transfer term q refers to the heat transferred to the fluid 

element from its surrounding fluid elements and that it is not equal to zero for the fluid element in the 
viscous region, even when the flow of the fluid as a whole is adiabatic, i.e., there is no heat transfer 
between the fluid and the bounding wall. For instance, in the preceding analysis of a two-dimensional 
boundary layer flow, it is seen that (i) there is a positive heat transfer into the fluid element in the 
boundary layer, and (ii) when the Prandtl number of the fluid is equal to one, this amount of heat 
transfer into the fluid element is equal to just the work done by the fluid element against the viscous 
forces acting on the fluid element, thereby keeping the stagnation enthalpy constant along the streamline. 
This is a very useful result, which provides a Bound basis in taking stagnation enthalpy or stagnation 
rothalpy constant along the streamline in current engineering calculation for turbomachine flows. 

Because of this canceling effect, it should be kept in mind that when the viscous effect on the stagnation 
enthalpy is considered in the calculation, it is not correct to keep one viscous work term in equation 
(2.20) or equation (2.53), or two viscous work terms in equation (2.23) and to neglect the heat transfer 
term in these equations. 

3. Entropy equation (2.10) or (2.61) clearly shows that the increase of entropy of the fluid element 
along the streamline is made of two parts, namely the heat transfer to the fluid element q and the work 
done by the fluid element against the external viscous forces as seen by an observer moving with the fluid 
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element, the dissipation function $ (see eq. (2.27)). Depending on the nature of problem, q may be 
positive or negative (usually positive), but $ is always positive. 

4. The heat transfer term q in the energy equation and entropy equation is heat transfer to the 
fluid element from its surroundings due to temperature difference. The viscous term in the energy 
equation and entropy equation is work done by the fluid element against the viscous force acting on the 
fluid element. They are independently evaluated according to their own definition, which is set up to the 
heat transfer term and work done term in the First Law of Thermodynamics (eq. (2.15)). In light of the 
preceding argument, it is easy to see that the frequently-appearing saying “work done by the frictional 
forces acting on the fluid element turns into frictional heat and is added to the fluid element,” 

dH/dt = F, • V, T ds/dt = F, V (for instance, eq. (3.19) on p. 51 of ref. 45 and equations in the 
middle of p. 277 of ref. 46) tire incorrect. 
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CHAPTER 3 


GOVERNING EQUATIONS FOR FLUID FLOW ALONG RELATIVE STREAM FILAMENTS 

3.1 Following Fluid Flow on Relative Stream Surfaces 

In order to solve the steady three-dimensional irrotational or rotational flow in a relatively simple 
manner, an approach was taken in reference 3 to obtain the three-dimensional solution by an appropriate 
combination of mathematically two-dimensional flows on two different kinds of relative stream surfaces 
(figs. 3.1 to 3.3). The first kind of relative stream surface is one whose intersection with a z-plane, either 
upstream of the blade row or somewhere in the blade row, forms a circular arc (fig. 3.1). The second 
kind of relative stream surface is one whose intersection with a z-plane, either upstream of the blade row 
or somewhere inside the blade row, forms a radial line (fig. 3.2). These two kinds of relative Btream 
surfaces were designated as stream surface Sj and S 2 , respectively. 

STREAM SURFACE OF THE FIRST KIND — Sj 

Shown in figure 3.1 is a stream surface of the first kind formed by fluid particles lying on a circular 
arc ab of radius oa upstream from the blade row. It is a generalization for three-dimensional flow 
from the cylindrical surfaced usually considered in the two-dimensional design of turbomachines. 

STREAM SURFACE OF THE SECOND KIND— S 2 

A stream surface of the second kind is shown in figure 3.2. The important surface of this family is 
the one that lies about midway between two adjacent blades, and divides the mass flow in the channel 
formed by the two blades into approximately two equal parts. This surface is designated as the mean 
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stream surface or mid-channel stream surface S 2m . For blades with all radial elements, such as the one 
shown in figure 3.2, it is convenient to consider a mid-channel stream surface formed by fluid particles 
originally lying on a radial line ab upstream from the blade row. Otherwise the radial line is chosen 
about midway in the passage with the fluid particles originally starting out from a curved line upstream 
from the blade row as shown in figure 3.3. 

In general, both of these two kinds of stream surfaces are employed in the solution of the three- 
dimensional problem flow field in turbomachines. The correct solution of one surface requires some data 
from the other, and consequently, successive solutions between the solution of one of these two surfaces 
are involved. Yet, the solution of flow on each surface is manageable with the efficient technique for 
mathematically two-dimensional problems. 

Relations Among Relative Velocity of Fluid, Coordinates 
of Stream Surface, and Normal to Stream Surface 


In general, the coordinates of the relative stream surfaces, the components of the unit normal n 
(figs. 3.1 and 3.2), and the velocity components are related by the following equations; 


s (r, <p, z) = 0 

(3.1) 

n r dr + n^ d^ + n 8 dz = 0 

(3.2) 

n r W r + n * W * + = 0 

(3.3) 
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3.2 Equation Governing Fluid Flow Along S 2 
Stream Surface of S 2 Stream Filament 

Because S 2 stream surface involved in the three-dimensional flow calculation is always a general 
twisted surface, whereas the Sj^ stream surface involved can be a surface of revolution, the S 2 stream 
surface will be considered first in the following treatment. 

When the fluid motion on S 2 stream surface is followed, equations (3.1) and (3.2) are used to 
eliminate one of the three independent variables, the <p coordinate. That is, any quantity q on S 2 is 
now considered as 


q = f [r, z, <p (r, z)] ( 3 - 3a ) 

The change in the quantity q along S 2 due to a small change in r while z is held constant is (see 
fig. 3.4) 


/ 



{dr 


dq d<p 
d<p dr ) 


dr 


(3.4) 


Substituting d<p/d r from equation (3.2), for dz = 0, into the preceding equation gives 

dq ^ dq _ n r 1 flq 
dr dr n^ r d<p 


(3.4a) 


in which the bold partial derivative sign is used to indicate this differentiation following the stream 
surface. 
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Similarly, 


= ^ i (3.5) 

dz dz r d<p 


Along streamline on S 2 

= w. + w. ^5 (3.6) 

Dt dr dz 

Continuity Equation 

When the fluid motion is followed along the S 2 stream surface and equation (3.4) and (3.5) are 
used, the continuity equation for steady relative motion becomes 


lf (^ + 3K) (3.7) 

r dr dz 


where 


C(r, z) = - + 




a W r 
dr 


+ n u 


dw v 

dr 


+ n. 


aw, 
dr ) 


This continuity equation is put into the following form: 


(3.8) 


d(rBpW r J + d(tB»W,) = o 
dr dz 


(3.9) 
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by using an integrating factor B, which is related to C by the following equation. 


D In B _ w din B + w dlnB = _ c (3.10) 

Dt T dr Z dz 


or 



(3.11) 


Equation (3.9) is the necessary and sufficient condition that a stream function ¥ exists. $ is 
related to velocity components by 


_ = rBpW, 

d r 


(3.12a) 


fi - rBpW 
dz 


r 


(3.12b) 


The difference in ¥ at two points j and k on the S 2 surface is 

*k _ = J*d¥ = J^rBp(W x dr - W r dz) ( 3 - 13 ) 

The preceding equation indicates that B is proportional to the angular thickness of a thin stream 
filament whose mid-surface is the stream surface S 2 considered herein and whose circumferential 
thickness is equal to rB. Indeed, if the mass flow into and out of the element of such a stream sheet (cut 
between two planes normal to the z-axis, and a distance dz apart (see fig. 3.5)) is equated to zero, and 
the distances dr and dz approach zero as a limit, the following equation is obtained: 


I II 
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(3.14) 


4'Wr) + 4^.) = 0 

dr dz 

Comparing this equation with equation (3.9) and considering the mass flow relations show r to be 
proportional to rB. This proportionality means that physically B is a quantity which is proportional 
to the angular thickness of a stream filament whose mid-surface is the S 2 surface considered herein. 

With this interpretation, B is immediately seen to be closely related to the angular distance between two 
adjacent blades. In actual calculation, only the ratio rB to (rB)j or r to r. is important. In general 
it is easier to obtain the variation in rB from the distance between adjacent streamlines obtained on 
surface than to evaluate B/Bj by equations (3.11) and (3.8) using data obtained on S 1 surfaces. 

It is seen from the preceding section that in following the fluid flow along a stream surface, a 
consideration of conservation of matter automatically changes the fluid flow on the stream surface to the 
fluid flow along the stream filament. In general three-dimensional flow, the thickness of the stream 
filament changes with respect to the two coordinates r and z. In the case of S 2 stream filament (see 
fig. 3.6) it is easy to see that: (1) in the radial direction, r increases with the radius as the 
circumferential distance between the adjacent blades increases in the radial direction and (2) in the flow 
direction, due to the blade thickness, t decreases as the fluid enters the blade channel and then increases 
as the fluid moves toward the trailing edge of the blade. 

Dynamic Equation 

For general rotational flow, the dynamic equation (2.69) in the three perpendicular directions are 


CD 

<fi> 
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dW T 

+ w. 

aw r aw,' 

= _dl + T de 

(3.15a) 

l dr 
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, dz dr t 
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(3.15b) 
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dW, 
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di T ds 

dz dz 


(3.15c) 


In following the motion on S, equntion (3.15) me reduced to the following forme by using equntione 
(3.3) to (3.6): 


w » a M + w - f^l 
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where F is a vector having the unit of force per unit mass of gas defined by 


F 


1 1 1 dp 
n^ p r dtf> 


n 


(3.17) 
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3,3 Principal Equation for Fluid Flow Along S 2 Stream Filament 

By considering the fluid flow along an S 2 stream filament and using the partial derivative of fluid 
quantities as functions of two independent variables r and z, the principle of conservation of matter 
leads to the continuity equation given by either equation (3.9) or (3.14), which are the necessary and 
sufficient conditions for a stream function 9 to exist. The relative velocity components Wr and Wz 
are related to the partial derivatives of ¥ by equation (3.12). When this relation is used, the dynamic 
equation in the radial direction can be used to form a principal equation governing the fluid flow on S 2 
stream filament. The solution of the fluid flow on S 2 stream filament is concentrated on solving this 
principal equation for the single dependent variable 9 . It is much better than solving a number of 
dependent variables from a set of partial differential equations. The form of the principal equation will 
be given for the direct problem and the inverse problem in the next two sections. 

Principal Equation for Direct Problem 

In the direct problem the shape of the S 2 stream surface is given. In practice the shape is specified 
by a number of coordinate points lying on the surface. The components of the unit normal are then 
calculated by equation (3.2) as follows: 

Along the intersection of the S 2 surface and a constant- z plane 

(3.18) 

and along the intersection of the S 2 surface and a constant r-surface 
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(3.19) 



Component, n,, n p , nnd n, we determined by the preceding two equntion. Jong with the following 
equation: 


n l + % + 1 


(3.20) 


Let 


n r Fr 



(3.21a) 


n » _ Fa (3.21b) 

n P ^ v 


From equation (3.3) 


w v = -(/iW r + I/W f ) 


(3.22) 


Using equations (3.21a) and (3.21b) and other basic relations, the dynamic equation in the radial 
direction can be transformed into the following form (ref. 48). 
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where 
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It is seen from the coefficients of the second order partial derivatives that when W > a or < a, 
equation (3.23) is hyperbolic or elliptic. 

Procedure of Solution for Direct Problem 

In the direct problem of fluid flow along the S 2 stream filament the variation of filament thickness 
t (or rB) relative to its inlet value and equation (3.21) are given. There are nine independent equations 
governing the fluid flow, namely equations (3.12a), (3.12b), (3.22), (3.16b), (3.21a), (3.21b), (2.67), 

(2.71), or (2.72). The nine independent variables to be determined are W r , W^, W g , F r , F^, F z , I, S, 
and «. (It may be noted here that, when equation (2.67) instead of the complete viscous equation (2.66) 
is used for the energy equation in the calculation, there are only three independent equations among the 
dynamic equations in three directions and the energy equation, because the latter can be obtained from 
the former and the normality conditions between F and W.) The procedure of calculation is as follows: 

(1) Starting from an estimated * field at the beginning of calculation or from the * field 
determined in the previous cycle, compute Wr and Wz from equations (3.12a) and (3.12b), 
respectively. 

(2) Compute from equation (3.22). 

(3) Compute + wr and then F^ from the dynamic equation in the circumferential 

direction, equation (3.16b). 

(4) Compute Fr and Fz from equations (3.21a) and (3.21b). 

(5) When the approximate equation, equation (2.67) is used for the whole flow region, stagnation 

rothalpy I is taken to be constant along all streamlines. 

(6) For an invicid isentropic calculation, entropy s is taken to be constant over the whole region. 
For an analysis and design calculation, which tries to approximate the real flow as closely as 
possible, a certain empirical variation along the streamline with the difference between the exit 
and inlet value given by equations (2.71) or (2.72) is considered. In transonic turbomachines the 
abrupt entropy increase across the shock is also taken into consideration. 

(7) Solve V from the principal equation (3.23). 
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Repeat calculations 1 to 7 until the desired accuracy is reached. 

Principal Equation for Inverse Problem 

In the inverse or design problem of fluid flow along the mid-channel S 2 filament the variation of the 
filament thickness r (or rB) relative to its inlet value is empirically determined by the desirable blade 
thickness distribution. (In three-dimensional solution the Sj filament thickness is taken from the 
solution of S l filaments obtained in a previous cycle.) Now there are only seven independent equations 
governing the fluid flow, namely equations (3.12a) and (3.12b), (3.16a) to (3.16c), (c), (2.67), and (2.71) 
or (2.72), two less than that in the direct problem. On the other hand, however, there are nine 
independent variables to be determined. 

The differential of the coordinates of the S 2 stream surface are related to the F components by 

F r dr + F^rdp + F.ds = 0 (3-24) 


In order for this differential equation to lead to an integral surface of the form represented by 
equation (3.1), F must satisfy the following condition of integrability: 

F • V x F = 0 ( 3 - 25 ) 


Writing equation (3.25) in scalar form and using the relations (3.12a) and (3.12b) gives 



(3.26) 



By integrating along a constant r-line equation (3.26) gives 



If Fr = 0 at z 0 



(3.27) 


(3.28) 


Thus, there is only one degree of freedom left to the designer. Of all of the appropriate ways of 
utilizing this degree of freedom, the one found most useful is to prescribe an appropriate variation of V 0 
or V^r on the S 2 m surface, i.e., the following equation is prescribed: 

V$r = G(r, z) (3 ' 29) 


The principal equation formed by combining the continuity equation and the dynamic equation in 
the radial direction is 
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(3.30) 
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From the coefficients of this second order partial derivative, the principal equation (3,27) is seen to be 

I 2 2 

hyperbolic or elliptic then the meridional velocity Wj = J W r + W f is greater or less than the speed 
of sound. 

In turbomachine design the S, m surface may be specified in some manner other than 
equation (3.29). For instance in the case of an axial-flow turbine, in order to design the cooled rotor 54 
blades with no or minimum radial twist, the following function may be specified on the mid-channel S, 
stream surface: 
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(3.31) 


^ = 8i(*) 

w. 


In the case of a centrifugal compressor, in order to design the impeller blades with minimum deformation 
and stress during operation, the mid-channel S 2 stream surface may be specified to consist of all radial 

elements, i.e M 

^ = rg,M < 3 - 32 ) 

W, 


In general then, the S 2 stream surface may be specified by the following relation: 


W p 

W. 


g (*> a) 


(3.33) 


By using equation (3.33) and other basic relations the dynamic equation in the radial direction can 
be transformed into the following form (ref. 48): 
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where 


48 


l 1 1 



M = _ainB + ^ 

dz dz 


dl 

dz 


W2 a, 

* 

a 2 a* 


N = 


- (l + g 2 ) 



/ ' 


r 1 

1 dlnB _ da* 1 

H + J,- wjg£l 

, a 2 - W 3 
+ g 

g + <|g 

r dr dr a * 

l®r * dr) 

a 2 

J dr, 


a 2 - W 2 


a 2 W 2 


+ T — + F r + 2wW w 
(dr dr * 


Procedure of Solution for Inverse Problem 

The procedure of calculation of S 2 m stream filament is as follows: 

(1) Starting from an estimated 9 field at the beginning of calculation or from the 9 field 
determined in the previous cycle, compute Wr and Wz from equations (3.12a) and (3.12b), 
respectively. 

(2) When (V^r) is prescribed, compute from V^. When (W^/Wz) is prescribed, compute 
Wp from equation (3.30) and Wz is computed from step 1. 

(3) Compute and Fz from, respectively, equations (3.16b) and (3.16c). 

(4) When the approximate equation (2.67) is used for the whole flow region, the stagnation rothalpy 
I is taken to be constant along all streamlines. 

(5) For an invicid isentropic calculation, entropy s is taken to be constant over the whole region. 
For an analysis and design calculation, which tries to approximate the real flow as closely as 
possible, a certain empirical variation along the streamline with the difference between exit and 
inlet value given by equation (2.71) or (2.77) is considered. In transonic turbomachines the 
abrupt entropy increase across the shock is also taken into consideration. 
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(6) Solve $ from the principal equation (3.30) or (3.34). 
Repeat steps 1 to 6 until the desired accuracy is reached. 


3.4 Principal Equations for Fluid Flow Along Stream Filament 

The S, stream surface, new the hub and c«ing wall, are usudly considered to be surfaces of 
revolution. If the radius of the hub wall increase, (or decreases) in the How direction while the radio, of 
the casing wall decrease, (or incr«a«.) in the flow direction, the 8, surface in the mid-span region may 
be .lose to surface of revolution. Just like in the case of S, flow, the continuity equation requires that 

the mass of fluid flow in a thin stream filament be conserved. 

The flow equation, expressed by a set of orthogonal curvilinem coordinates (., p) (H- »■•) « " 
follows (refs. 2 and 4): 

Continuity Equation: 


a(rpW f ) + a(rpW y ) = 0 

c« dip 


(3.35) 


$ ,— W relations 


£ = 
ct 


(3.36a) 
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(3.36b) 


= vw, 

r dp 


Dynamic Equation: 
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(3.37) 


Principal Equation: 
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For general Sj surface or filament, where the surface twists in the circumferential direction, one 
may formulate the governing equation with respect to the coordinates as follows: 

(1) For turbomachines with axial inflow and axial outflow any quantity q on the S 1 surface is 
considered a function of <p and z. 

(2) For turbomachines with radial inlet and radial discharge, any quantity on the Sj surface is 
considered a function of p and r. 
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(3) For turbomachines with axial (or radial) inflow and radial (or axial) outflow, any quantity on 
the S 1 surface is considered a function of (<p, t). 

(4) Same as that in case (3) but considered a function of noncurvilinear coordinates x 1 and x 2 . 
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CHAPTER 4 


THREE-DIMENSIONAL FLOW EQUATIONS EXPRESSED IN TERMS OF 
NONORTHOGONAL CURVILINEAR COORDINATES AND CORRESPONDING 
NONORTHOGONAL VELOCITY COMPONENTS 

Equations governing the fluid flow along Sj and S 2 stream filaments were first derived in 
references 2 to 4 employing orthogonal curvilinear coordinates (r, <p, z) or (p,l, n). These equations have 
been adopted in treatises on turbomachinery (for instance refs. 45 to 47), programmed into computer 
codes (for instance, refs. 8, 10 to 19, 22, and 23) and used in analysis and design of turbomachines (for 
instance, refs. 21 to 26). During calculation of actual engineering problems it soon becomes evident that 
in order to improve the relatively low accuracy of numerical differentiation (ref. 27) occurring at grid 
points, unequally spaced near a curved boundary the coordinate line may coincide with the leading and 
trailing edges of the blade (see fig. 4.1) and a computer code may be used universally for turbomachines 
of different geometry. Work to employ general nonorthogonal curvilinear coordinates for the calculation 
of fluid flow along Sj^ and S 2 stream filaments began in the sixties (refs. 28 and 29). Presented in the 
following sections are some basic relations of general nonorthogonal curvilinear coordinates, equations 
governing fluid flow along Sj and S 2 stream filaments employing general nonorthogonal curvilinear 
coordinates and the methods of solution. 

4.1 General Curvilinear Coordinates in a Three-Dimensional Space 

Let (y 1 , y 2 , y 3 ) be the usual orthogonal Cartesian coordinates of point P and (x 1 , x 2 , x 3 ) be its 
nonorthogonal curvilinear coordinates. In general the square of the element of arc ds is 


53 



(4.1) 


ds 2 = gy dx 1 dx* (i, j = 1, 2, 3) 


where gy is the covariant metric tensor of the three-dimensional space 


gij = £L £L = Sj, (k = X, 2. 3) 

dx 1 dx? 


(4.2) 


Sij = fa cos *ij 


(4.3) 


where fly is the angle between the two basic vectors «j and «j. 

The lengths of the elements of arc measured along the coordinates lines of our curvilinear systems, 
are 



(4.4) 


The element of volume dv is 


dv 



dx 1 dx 2 dx 3 


(4.5) 


where g is the determinant |gy|. Let e‘ denote the reciprocal base vectors defined by 
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( 4 . 6 ) 


«1 ' 


J = «< or e 1 = 1 X * 

£ 


The corresponding contravariant metric tensor of the three-dimensional space is 


g u = ^ 

g 


= e 1 • e* = g* 1 


( 4 . 7 ) 


where is the cofactor of the element gy in g. 


For both systems it is convenient to use base vectors of unit length defined by 


«i = ®l * 



( 4 . 8 ) 




( 4 . 9 ) 


A vector B in the three-dimensional space is now either expressed in terms of the base vectors €j 
and U|, or the reciprocal vectors e* and u 1 , as follows: 



( 4 . 10 ) 
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B = B'uj 


( 4 . 11 ) 


B = bje 1 


( 4 . 12 ) 


B = Bju' 


( 4 . 13 ) 


where b' and B’ are, respectively, the contravariant components and the physical components along 

of the vector B; and b; and B; are, respectively, the covariant components and the physical 

components along e 1 of the vector B. The covariant component and its corresponding physical 
component can be calculated from the contravariant component with the following formula: 


bj = B • «j = gy b* 


( 4 . 14 ) 


Finally the differential operators in general curvilinear coordinates are as follows: the gradient of a 
scaler I is given by 


VI 



( 4 . 15 ) 


the divergence of vector W is 
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V . W = 1 a ^ wi ) 

f* 


(4.16) 


and the contravariant components of vector E = V x V = £*ej are 
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dv 3 dv 2 
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k dx 2 dx 3 
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dvi dv 3 
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(4.17) 
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4.2 General Curvilinear Coordinates on a Surface 


In the investigation of fluid flow along an Sj or S 2 stream filament the governing equations are 
written for the fluid flow on the mid-surface of the filament. Such a surface is a two-dimensional 
manifold embedded in a three-dimensional enveloping space and is usually described by coordinates u a 
(a = 1, 2), called the curvilinear or Gaussian coordinates on the surface. 

Under this notation the relevant equations corresponding to those in the preceding section are as 
follows: 
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( 4 . 18 ) 


ds 2 = & a p dx° dx^ (a, $ = 1, 2) 


where a a ^ 


is the covariant metric tensor of the surface 


dy k ay k 

a = IfL JL = e a • e. 


= a 


du° 


Pa 


( 4 . 19 ) 
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aP 



cos e afi 


( 4 . 20 ) 


The lengths of elements of arc measured along the coordinate lines are 



( 4 . 21 ) 


An element of area dA is 


dA = 



( 4 . 22 ) 


where a is the determinant |a a ^| . 
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= \* a p\ 


— a 


11 a 22 


sin 2 0 12 


( 4 . 23 ) 


I I! 
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The contravariant metric tensor of the two-dimensional surface is 


a°^ = e° • J = aP a 

a 11 = a 2 2 -r a, = (a n sin 2 0 12 ) 1 ( 4 - 24 ) 

a 22 = a n -f a = (a n sin 2 (? 12 ) 1 


.12 


.21 


*12 


a = 


cos 


^12/ 


1 \/ a ll a 22 


sin 


'12 


(4.25) 


In the investigation of Sj and S 2 flow employing nonorthogonal curvilinear coordinates the two 
Gaussian coordinates, u 1 and u 2 , on the surface, are selected to be the same as two of the three 
curvilinear coordinates in the three-dimensional space. For instance, in the case of Sj flow, x and x 
shown in figure 4.1 are two nonorthogonal curvilinear coordinates of the three general curvilinear 
coordinates referring to the three-dimensional space and are the same as the two Gaussian coordinates 
which refer to the two-dimensional surface. The third coordinate x 3 , which refers to the three- 
dimensional space, may be selected as normal to the x 1 — x 2 surface. In that case W 3 = w 3 

= W 3 = w 3 = 0 - 

In the following, distinctions between (u 1 , u 2 ) and (x 1 , x 2 ) will not be made, but e ft , a . . . (a, 
/? = 1, 2) will be used for the Sj surface, whereas e*, gy, g . . . (i, j = 1, 2, 3) will be used for the three- 
dimensional space in which the surface is embedded. For a given problem, are computed from 
equation (4.19) through numerical differentiation. 
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For vector W on the S t surface the following are obtained: 


W j = W ■ ej = a u w 1 + a 12 w 2 = ^(w 1 + W 2 cos 0 12 ) . 
w 2 = W • e 2 = a 21 w 1 + a 22 w 2 = f^( W, cos 9 n + W 2 ) 


(4.26) 


W! = ^"w! = (w 1 + W 2 cos 0 12 ) -T sin 0 12 
W 2 = = (W 1 cos 9 n + W 2 ) -h sin 9 U 


(4.27) 


The preceding equations can also be obtained by the geometrical relations indicated in figure 4.2. 

4.3 Basic Equations Governing Fluid Flow on Sj Surface of Revolution 

General Equations: Denoting the normal distance between two adjacent surface of revolution by r, 
the continuity equation of steady relative flow along the stream filament is 


V • (rpW) = 0 


(4.28) 
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Substituting equation (4.16) into (4.28) results in 


djrpW 1 J&22 8 in ^12) d\rpVt 2 sin 0 12j _ ^ 
ax 1 


(4.29) 


ax - 1 


This equation can also be obtained by considering the mass flow into and out of the elementary volume 
(ds 1? ds 2 , r) (see fig. 4.3). 

When expressions (4.10), (4.15), and (4.17) are used, the dynamic equation (2.53) in the e 2 
direction is 


-w^! x ?e 3 = 


di_ _ T jte 
{ ax 2 dx 2 


(4.30) 


When equation (4.10) and the following relation between the absolute and relative vorticities are used: 


V x V = V xW + 2w 


(4.31) 


the dynamic equation in the e 2 direction for steady relative flow is (fig. 4.4). 
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(4.32) 


-w 1 

( ^2 _ 

dwj' 

- 2wfi sin a 


_ T ds 


.U* 1 

dx 2 ; 


dx 2 

dx 2 


By using equation (4*26) the preceding equation becomes 



This equation can also be obtained by applying Stokes theorem to the surface element of figure 4.3. 
Using relations (4.10) and (4.15) the energy equation (2.51) for steady relative flow becomes 


5! = W • VI = w^ • 
Dt 


« e l = w> il = 0 

dx 1 dx 1 


or 


m _ _ai_ + wl dI = n 

Dt fin dx * fi*i dx * 


Special Forms of Flow Equations on S^: 


(4.33) 


(4.34) 


(4.35) 
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A number of useful special forms of equations presented in the preceding section can be obtained by 
selecting a number of different x 1 - and x 2 -coordinate lines, as follows: 

(a) x 1 -lines coincident with streamlines. In this case W 2 = w 2 = 0, W 1 = W, equations (4.29), 
(4.33), and (4.35) simplify to, respectively, 



(4.29a) 


d 

d ^ 2 



d 

dx 


_ (w 1 \j &22 cos ^ 12 ) ~ 2 ^a n a 12 ta sin a sin 0 12 


+ 



_ai_ _ T _a^ 

,dx 2 dx 2 , 


(4.33a) 


(4.35a) 

If at the same time x 2 is taken to be normal to x 1 , the preceding equations are further simplified with 
cos 0 12 = 0 and sin 0 12 = 1. 

(b) xMine coincident with the meridional coordinate line ( and x 2 = <p. In this case, 0 12 = 270°, 

w 1 = W„ W 2 = - W p , = 1 , = r, equations (4.29), (4.33), and (4.35) become those used in 

references 4 to 9. 

(c) Sj surface is a cylindrical surface of radius r. If y = x 2 = Tip and x 1 = z, equations (4.29), 
(4.33), and (4.35) become those of plane flow as used in reference 54. If, however, xMine is taken to be 
coincident with the streamline W 1 = W, W 2 = 0 and the equations are further simplified to 


31 

dx 1 


= 0 
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(4.29c) 




(4.33c) 


_£! = 0 (4.35c) 

dx 1 

Furthermore, if x 2 -line is taken to be normal to x 1 -line, an interesting result is obtained the quantity 
W is an invariant along an x 2 - line in the case where the variations of I and s along the x - 

line are equal to zero. 

As seen in the preceding equations, dynamic equations (4.33) and (4.33a) contain the variations of I 
and s of the fluid in the e, direction. When they are uniform at the blade inlet and no change occurs 
along the streamline, they are uniform everywhere— this is usually assumed on S x flow calculation and 
will also be assumed in the following discussion of methods of solution. In cases where their derivatives 
are not very small, they should, of course, be properly considered. 

4.4 Methods of Solution for Sj 

When the x 1 -lines are corrected during calculation to coincide with the streamlines, equation (4.29), 
along the xMine, gives 
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(4.36) 


pW = G>W)j 



sin 0 


A 


7 



sin 0 12 


in which ^a 22 and 0 12 are computed after each correction of the streamline coordinates. For 
nonisentropic flow of a perfect gas 


p _ 

r \ 

h 

1 8-Si 

1 e = 

, _ 7 - 1 0>W) 2 

Pi 9 ' 

H? 


2 (/»i tt i) 2 (/»//»i) 2 


(4.37) 


The values of p and W can be determined from their product by rewriting equation (4.37) as (ref. 3) 


£ 2 = 


/ \ 
i - L 
l £ J J 


2 

T 3 ! 


(4.38) 


where 


E= L e R 


(4.39) 


65 




p w 




- 2 . 


8 — 8 : 


R 


The table in reference 3 can then be used find E (i.e., p) from the known value of $ (i.e., pW). W is 
then computed from pW after p obtained. 

When the derivatives of I and s are neglected, equation (4.33a) yields 

— ( w l/aii1 = — ( W f 22 cos ^12) - 2 ^ a u a 22 w 8in * sin *12 ( 4 - 41 ) 

ax 2 dx 1 


This equation is used to obtain a new set of values of W at the grid points by integrating along the 
x 2 -lines with the right side of the equation computed from known values of the preceding cycle. The 
integration may be started from the mid-channel streamline as in references 4 to 6 or from the stagnation 
streamline. The mass flow along the x 2 -line is then computed by 

m = f* 2 rpW sin 0 l2 1 /a^dx 2 ( 4 ‘ 42 

A comparison between the computed mass flow and the correct mass flow leads to new values of 
streamline positions. This process of correction should be applied at each station along the streamlines, 
quite similar to the method employed in references 4 and 6. It may be noticed, however, that with the 
use of the present nonorthogonal curvilinear coordinates the tangency of the relative fluid velocity at the 
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blade surfaces is automatically satisfied. This method, of course, can also be employed to modify existing 
blading or to design new blading as in references 4 and 35. 


If it is desirable to keep the xMine fixed during calculation, then 

(1) pW 1 along the xMine is corrected by the complete continuity equation (4.29). 

(2) compute />W by 


pW 


(pW 1 ) 2 + t>W 2 ) 2 - 2kw 1 )G>W 2 )cos 



Vi 


( 3 ) 


integrate (w 2 ^a u j along x 1 -line by 


i( wl f^) _ 


3x 2 


= — 2 ^a u a 12 u sin a sin 


'12 





(4) compute mass flow along 


x 2 -line by 



sin 0 


12 \ a 22 


dx 2 


(5) correct the W 1 values at the grid points in the same manner as in the previous case. 


(4.43) 


(4.44) 


(4.45) 
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The method described above is quite similar to that used in references 4, 6, and 54 except that the 
solution of the latter is extended in the circumferential direction from the mean streamline by using a 
Taylor series consisting of three terms. Whereas, the solution of the former is extended circumferentially 
to successive streamlines by a Taylor series of only two terms. 

In general the continuity equation and the dynamic equation are to be solved simultaneously through 
the use of the stream function ¥. Equation (4.29) is necessary and sufficient condition that a function $ 

exists with 


if. = - r/>W 2 ^7 s » n ^12 


dx 1 


if. = rpW 1 sin 0 l2 

o * 


dx 


(4.46) 


Substituting equations (4.46) into (4.33), the following principal equation governing the fluid flow along 

X 2 • 

the S t filament of revolution expressed in terms of nonorthogonal curvilinear coordinates x - and x is 

obtained. 


i d 2 * _ 2 cos *12 a 2 * + J_ a 2 * + _J_ a* + Jf_ a* _ M ( 4 . 47 ) 

a n a^i ) 2 fiuH*** 1 *** a22 d(x 2 ) 2 fin 0x1 fin 9x2 


where 
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I'll: 




K = - 



cos 7 din r 

+ + 

^a n dx 1 



d0i2 

dx 1 


r 

1 dln/»_ 

1 

cos 0 12 din p 

i a* + 

1 dlnp cos ^12 dlnp 


/ a ii~ ^ x , 


^ dx 2 ^7 3x1 , 


-I- * ^ — 2wr/» sin a sin 2 0 12 

^2 


Equation (4.47) may be considered as a generalization of the principal equation (3.38) in which the 
orthogonal curvilinear coordinate. (I, p) me u«d. When x 1 and X 1 become 1 and p, respectively, 
equation (4.47) reduces to equation (3.38). 

When the ¥ -derivatives are evaluated with suitable numerical differentiation formulas, including 
those at the unequally spaced grid points (ref. 27) the resulting set of a large number of algebraic 
equations involving the unknown ¥’s at the grid points may be solved either by the direct matrix 
method (refs. 3, 4, 8, 10 to 12, and 27) or by the relaxation method (refs. 3, 4, 10, 11, and 27). The 
boundary conditions of the periodic variation of flow and the flow angles at the inlet and outlet stations 
are to be satisfied in the same manner as described in references 4 and 8. Inasmuch as fixed (x 1 , x 2 ) 
coordinates are generally used in solving the principal equation (4.47), J and K are then a function 


69 



fixed geometry and have to be computed only once. The exclusion of derivatives W 1 and W 2 on the 

right side of equation (4.47) is helpful in the convergence of the solution. 

4.5 Equations Governing Fluid Flow on a General S 2 Surface 


Vorticity Equation and Dynamic Equation: In general, the problem of fluid flow along S 2 stream 
filament is treated in the same manner as described for S x stream filament m the preceding sections. 

But there is an important difference in the S 2 surface inverse problem in that the shape of the S 2 
surface is known only after the solution is obtained. Because of this and also for the elimination of the 
angular derivatives in the governing equations, the independent variables, i.e., the nonorthogonal 
curvilinear coordinates x 1 and x 2 are selected to lie on the meridional plane (see fig. 4.5) while 
following the fluid motion on the S 2 surface in exactly the same manner as was done m references 3, 10, 
and 11 in which orthogonal coordinates z and r on the meridional plane were used. The third 
coordinate x 3 for three-dimensional flow is chosen to be the angular coordinate <p, which is the same as 
that used in the cylindrical coordinate system previously employed in references 3, 10, and 11. 
Corresponding to this choice of x 3 , ^33 = r, g = r a, )[g — r , W — W^, w — . r, V 

V tf , v 3 = V $ + r, W 3 = W^, w 3 = W^r, V 3 = V 0 and v 3 = V tf r. 


It was found in references 2 and 3 that it is more convenient to use the absolute vorticity than to use 
the relative vorticity in treating the relative flow on S 2 surface. Under the present nonorthogonal 
curvilinear coordinate system, the three contravariant components of V x V (eq. (4.17)) becomes 



By using these expressions, the three nonorthogonal dynamic equations in the ej directions are 


w 2 

dw 2 dwj 

w 

dw x d(Vfr) 

dl _ T ds 

\j & 22 

dx 1 dx 2 

r 

dp dx 1 

dx 1 dx 1 


w 

d(v,r) _ dw 2 

_ W 1 

dw 2 

dwj 

dl _ T ds 

r 

. dx 2 d(p 

/ a n 

dx 1 

dx 2 . 

dx 2 dx 2 


W 1 

dwj 

a M 

W 2 

d(v) dw 2 

= dl _ T ds 


dip 

dx 1 

\[*22 

dx 2 

dip dp 


(4.48) 


(4.49a) 


(4.49b) 


(4.49c) 
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Alternatively these equations may also be obtained by using the following expressions of the rectangular 
components of W and V x V in the e\ e 2 , and e 3 directions (see fig. 4.2). 


W c = W l sin 6 U , W d = W 2 sin 0 12 , W„ = W p 


(4.50) 



(4.51) 


We notice that the form of equations (4.49) and (4.51) are quite similar to equations (19) and (75), 
respectively, of reference 3 employing orthogonal curvilinear coordinates (z, r, <p), but that the w’s are 
now more complicated covariant components of the relative velocity W. 

With different choices of x 1 - and x 2 -coordinate lines suitable for different problems, the general 
form of equations (4.49) reduces correspondingly to a number of special forms. For instance, if it is 

selected that 

(1) x 1 -line is coincident with the projection of the streamline on the meridional plane, then 

W 2 = 0, W 1 = Wj; and w x = s] & n Wj, w 2 = W J &22 COB ^12> an< ^ e 9 uat, i° n ( 4 - 4 9) 

becomes 
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w 

V 

s(w,^) 

r 



3 (v»r) 

di 


ai T ds 

d~t di 


w 

<p 

a(v) _ 

a(w. ^22 cos ^12) 

W, 

a(w ( ^cos # 12 ) a(w ( ^) 

ai _ T as 

r 

. dx 2 


17 

a x 2 

ax 2 ax 2 


(4.52) 


W| 




di 


£_ T as 

dip dp 


2 • 

(2) xMine is coincident with the projection of the streamline on the meridional plane and x -line to 
be radial line normal to the z-axis, then = a, 0 2 — 0, 0^ 2 = 90° — a, 
a = Ja n a 12 cos a, w 2 = Wj ^a 22 sin <r, equation (4.51) becomes 
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e = 


1 

d 

M 

d(W| ^7 sin <r) 


a u a 22 cos (7 


dr 

dp 


1 

s(w ,f^l _ a(v) 


r cos ° 

dp o< . 


| 

1 

a(w, /a22 sin a) 

a(w,^) 

r ^a ^ 

di 

dr J 


(4.53) 


and (4.49b) becomes 


s(v) _ 

d(w | ^22 s * n ^ 12 ) 

w« 

s(w, 

^ a 22 sin ^12) 

d(w| ^7) 

_ di _ T da 
dr dr 

( 4 . 54 ) 

dr 

dp 



o< 

dr 



and (4.49,) M d (4.49c) remain unchanged. Equation (4.54) may be used to evaluate the change of W, 

along the r-coordinate lines for axial-flow machines. 

Under this coordinate system, the correct expressions of V x V in the r, p and z directions in 

terms of the rectangular components Wj and or are 
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_ tan a 

1 a(v,r] 

l_ 

i 

a M . 


1 3W, 

r 


sin o 

dl 

sin <r d<j> 

S3 — 

tan a 

i 

i 

sT| 

1 

aw, 

~<t> 

£ : 31 

r m cos o y^a 2 2 

dr 

H W ^ 

1 

cos a = _ 

i 

a M _ 

d(w,sm 



Z “ 

r 


dr 

dp 




dWjsin « 

d 4 > 


equation (4.51) reduces to 


line and 

9 X = 

2 * 0, then 

Wj = W 1 

fin* 

A 1 = _ 

1 

a(v) _ a(w 2 


> 

/ a lla a2 

. 5x 2 

a* 


e = - 

1 

a(w‘^] 

_ a M 

► 

r i 

/ a ll a 22 

d<t> 

dx 1 


^ = _ 

1 

a(w 

_ a(w l 

r i 

/ a ll a 22 

dx 1 

dx 2 J 


and (4.49) reduces to 


(4.55) 


a ll a 22 > 


(4.56) 
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w 2 

a(w _ s(w l f^) 

w 

_ V 

a(w l fi^j _ a(v) 

ai _ T as 

fin 

dx 1 dx 2 

r 

d<P dx 1 

dx 1 dx 1 


w„ 

a(v^) _ a(w 

_ W 1 

s(w 

1 _ a(w'^) 

_ di _ T ds (4.57) 

r 

dx 2 

dip 

fin 

[ 

dx 1 

dx 2 

dx 2 dx 2 


W 1 

a(w‘^7) _ 

a(v) 

W 2 

a(Vtf r ) _ a(w 2 fi^ t 

di _ T ds 

fin 

dp 

dx 1 

fin 

dx 2 9< p 

dp dp 


If e n = 90° and 9 X = 0 2 = 0, then ^7 = fin = 1 > equations (4.55) and (4.56) reduce to, 
respectively, (19) and (95) of reference 3 in which the orthogonal cylindrical coordinates r, <j>, i are used. 

(d) x 2 -line is normal to x 1 - line and xMine is coincident with the projection of streamline on the 
meridional plane, then = 0 2 = W 2 = ~ w l ~ fill' w 2 ~ ^ ~~ filfiii ant ^ 

equation (4.48) simplifies to 


76 


1 





■/ 


r i/a n a- 


ll a 22 


dx 2 



(4.58) 


and correspondingly 


S 1 _ 

i sM 


r \^ a 22 

-2 1 
M 

a(w _ 3(v,r) 

fin 

dp c< 

5 3 _ _ 

i a(w'v^ 


l/ a ll a 22 


It is interesting to notice that in the present coordinate system the components of the absolute vorticity 
H 1 &nd contain only one term each, and that for flow, in which the absolute vorticity remains zero, 
both V^r and W 1 remain constant along the x 2 - line, while W 1 must vary with respect 

to the angular coordinate <p according to the variation of V with respect to the meridional 
coordinate I as the fluid moves through the blade passage. In this coordinate system, equations (4.49) 
simplify to 
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w 

a(w, f^) 

| _ 3| 

M 

_ 31 _ T 3s ( 4 - 60a ^ 

r 

[ 

dp 


3t 

31 31 


w„ a(v,r) _w^ a(w,^ = 31 _ T^l (4.60b) 

r 3x 2 ^a n 3x 2 3x 2 3x 


W| 

a(w‘ 

1 

QJ 

31 _ T 3s (4.60c) 



3^> 

3* 

d<p d<p 


a n 


4.6 Governing Equations for S 2 with Independent Variable <p Eliminated 


The unit vector normal to the S 2 surface is expressed as 


n = nje 1 = NjU 1 


(4.61) 


the orthogonal relation is expressed by 


N, 


No 


nJ ■ dx‘ej = n^x 1 + n 2 dx 2 + n^d ip = — Idx 1 + — L d * 2 + = 0 

.Lll 2 


(4.62) 
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The derivative of fluid quantity on S 2 stream surface with respect to x 1 or x 2 with the other kept 
constant is 


± = ± + ±i* = ± + ±**=±-?i± = ±- — l—± (4.63a) 

dx 1 dx 1 dx 3 dx 1 dx 1 d<P dx 1 dx 1 n *> d P dx 1 


d __ d + d <Jx 3 _ d d dtp __ d _ n 2 d __ d 
d? di 2 dx 3 dx 2 dx 2 dx 2 dx 2 n *> dx 2 


n 2 a 


V 



(4.63b) 


When equations (4.63) are used the continuity equation for steady relative flow becomes 



(4.64) 


where r is the circumferential distance between two adjacent S 2 surfaces and the dynamic equations 
become 


W 2 


dw 2 



dwi + d(Vgr) dl _ T da 
dx 2 r dx 1 dx 1 dx 1 


(4.65a) 
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(4.65b) 


vi, f a(v) _ 

W 1 | 

dw 2 

dw x ' 

dl _ 

Til 

r dx 2 


k dx 1 

dx 2 j 

d^ 

dx 2 


-h 


D (v) _ w 1 + Wl^M — f = F r (4.65c) 

Dt 9x1 ^*22 3x2 


Alternatively, these equations can be obtained by projecting the three components of equation (96) of 
reference 2 in e 1 , e 3 , e 3 directions. F in these equations is the same F of reference 2, fj being its 
covariant components 


F 


1 1 dp 
P d<p 


\ 

n 

j 


1 1 dp n = = FjU 1 

N^r p dip 


(4.66) 


From equations (4.26), 


£5 - = -£-[(w l cos e i2 + w 2 )fn\ - -^-[(w 1 + W 2 cos e i2 )f^] (4.67) 

dx 1 dx 2 dx'" dx 2 
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III 


Substituting equation (4.67) into equation (4.65) gives 



(4.68a) 


(4.68b) 


D(V) = djVjr) d(v,r) 

Dt fin 9x1 fin 9x2 




(4.68c) 


Equation (4.68b) is the most important equation for S 2 flow just as equation (4.33) for S t surface 

flow. 

Corresponding to the different choice of the x 1 - and x 2 -coordinates, equations (4.64) and (4.68) 
reduce to the following special equations: 

(1) xMine coincides with the projection of the streamline on the meridional plane. 

W 2 = 0 and equations (4.64) and (4.68) simplify to, respectively 
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JS * 

m 



(4.69) 


and 



W <p a 0V) _ + = (4.70a) 

~ JT~ a* ai 1 


_ 

1 

5 \ W I V a 22 008 *12 

l/^n 

w <t> Wf) ai 

- T * - f. 

dx 2 

dt 

w ( 

r dx 2 dx 2 

ax J 


(4.70b) 


D(Vjr) W ( d(V,r) _ _ 

= — — r ^> r 

“ »« 


Equations (4.69) and (4.70b) of S 2 surface are quite similar to equations (4.29a) and (4.33a), 
respectively, of Sj surface. 

(2) xMine coincides with the projection of the streamline on the meridional plane and x 2 -line is in 
the radial direction normal to the z-axis. Equations (4.69), (4.70a), and (4.70c) remain 
unchanged, but (4.70b) becomes 


PAGE IS 

*° 0R quality 
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(4.71) 



w, £(V) + ai _ T * 

r dt dr dt 


1 W 9 m 2 

(3) x-line coincides with the projection of the streamline on the meridional plane and x -line is 
normal to xMine. Equations (4.64), (4.70a) and (4.70c) remain unchanged, but (4.64) simplifies 
to 


a ( w i fin\ _ fin 

_ 3(V„r) + _W _ T as _ , 

0x 2 W, 

r 0X 2 dx 2 dx 2 2 


(4.72) 


It is noticed that equation (4.70b) contains an f 2 term which is, in general, not equal to zero inside the 
blade passage. 

By the use of equations (4.63), the energy equation (4.35) becomes 


DI W 1 dl W 2 01 = Q 

Dt fin 9x1 fin 9x2 


For the special case where xMine is coincident with the projection of the streamline on the meridional 
plane it simplfies to 


J?L = 0 (4.74) 

dx 1 
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4.7 Methods of Solution For S 2 Flow 

Since the editions governing the flow along S, stream filament are quite similar to thoee along S, 
stream filament (the main difference being the difference in the partial derivatives involved), the methods 
of solution are also quite similar. In the first group of methods of solution, W> or * is integrated along 
successive x’-coordinate lines, whereas in the second group of methods of solution «-values at all grid 

points on the mid-S 2 surface are solved at the same time. 

In the first group, if the x'-lines are corrected during calculation to coincide with the projection of 

the streamlines on the meridional plane, the continuity equation (4.64) gives 


pW, = (pW,), 



(4.75) 


S, m.—n. flow is calculated by considering the effect of entropy increase in the density equation 


s = JL 

o 


H - - V] 
2 9 


H: 


7— 1 8 — 


81 


* = 


1 - 1 
2 


\2 


pW| 

o o 

ITlM 


Hi 


Hi - i v; 


-r+i ,«-« 
-r-i ga R 


(4.76) 


The values of p and W can be determined from their product by rewriting equation (4.75) in the form 

of equation (4.38) with 


84 
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(4.77) 



a(V»r) + _ x _ f 2 

r dx 2 dx 2 dx 2 


The table in reference 3 can be used to find E (or p) from the known value of * (or pW ( ). W, is then 
computed from p Wj after p is obtained. 

The right side of equation (4.70b) is more complicated than equation (4.41) of S x flow. In the 
inverse or design problem, the variation of V^r is specified by the designer, whereas in the direct or 
analysis problem it is obtainable from the given shape of the S 2 surface (see next section for further 
discussion). The variation of I is obtained from its inlet value and equation (4.73). The variation of s 
is estimated from its inlet value and the empirical equation. W ( is then integrated along the x 2 - 

lines from a streamline along the inner boundary or in the middle of the flow region in the same manner 
as that in the S, surface flow calculation. In the case where the x 2 -line is chosen normal to the xMine 

it is convenient to rewrite equation (4.72) as 

m = KrffV rpWj sin 6 n f. ^ dx 2 ( 4 ' 7S 

u (* ) h 


and integrate (W, f^) 2 along the x 2 -lines. After variation of W, along x 2 is obtained, the mass 
flow along x 2 is computed. 


with 


K c = 


2xr - Nt 


2irr 

T 


/ ^ 


> 


/ \ 

i-22. 
1 2xr, 

_ 2irr 

i-l 

pJ 

_ 2irr 

p - 1 

~~ T 

~ T 

p J 


(4.79) 
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where { is a correction factor as indicated by figure 19 of reference 7. K the range of integration 

includes the wall boundary-layer regions, m calculated should be slightly greater than the actual mass 
flow because the calculated velocity in the wall boundary-layer region is greater than the actual value (see 
fig. 2.5). In the design problem we may either specify the mean streamline and obtain the coordinates of 
the inner and outer boundaries from given mass flow, or specify one boundary and obtain the other 
boundary. In the analysis problem, the mass flow calculation is used to correct the values of W, and p 

on the streamlines (see ref. 6). 

If it is preferred to use fixed xMines, then the variation of pW 1 along the x 2 -lines should be 
corrected according to continuity equation (4.64). Compute pW by 

\\ (4.80 

pW = {[(pW 1 ) 2 + (/>W 2 ) 2 - 2(/»W 1 )(pW 2 ) cos 0 12 ] + (/>W p ) 2 } 2 


integrate W 1 ^aj~j along x 2 -line by 


a( w1 ^) _ 
dx 2 w 1 

w p d(V e i) + di _ T ds_ _ 
r dx 2 dx 2 dx 2 


+ -t 

dx 1 

(w 1 cos e n + w 2 ) f^\ - — t [w 2 v 

dx 

£7" cos e l2 


(4.81) 


compute mass flow along x 2 -line by 


4f : 


.d-3 
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(4.82) 


in = K c 



sin 0 12 



and Correct W 1 values at the fixed grid points in the same manner as in the previous case. 

Analogous to the Sj-surface problem, the continuity equation and dynamic equation may be solved 
simultaneously through the use of stream function. From equation (4.64) a stream function can be 
defined with 


jj! = - rp W 2 ^a n sin 0 l2 
dx l 

— = rpW 1 f 2 7 sin 0 l2 
dx 2 


(4.83) 


Substituting equation (4.83) into the dynamic equation (4.65b) gives the following principal equation for 
flow along S 2 stream filament. 

± J?L - 2 _£?_ + ±_^L + _i_il + JLii = M (4.84) 

a ll d(x') 2 ^ x ** <a 632 ®( x2 ) 2 ^a u dx 1 ^21 9x 2 


where 
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M 


1 d In p _ cos e n d In p 1 d9 + r 8in 9 _w p Q 
J&22 ^ \[*U ^ ^ / a ll a 22 




Equation (4.84) may be considered as a generalization of the principal equation (107a) of reference 3 
where orthogonal curvilinear coordinates r, <p, z are used. When x 1 and x become z and r, 
respectively, equation (4.85) reduces to equation (107a) of reference 3. 

When the ^-derivatives are evaluated by suitable numerical differentiation formulas of equal or 
unequal spacing (ref. 27), the resulting set of a large number of algebraic equations containing the 
unknown 9-values at the grid points are obtained. Also the dynamic equation may be written as 


^t • WdL = J J (V x W)^ ^a n a 22 sin 0 12 dx l dx 2 


= f f ^7 K £(V) _ 

wl I r 0 X 2 Qx* 


(4.85) 


_5ii. + t __ — f 2 dxMx 2 


and integration along 1, 2, ..., 8 , 1 around the grid point E results in an algebraic equation containing 
the 9 -values at the grid points A,B,...,I, (fig. 4.7). The resulting set of algebraic equations containing 
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unknown ^-values at the grid points in the order of several hundred for one stream filament, can be 
solved by either the direct matrix method or relaxation method as employed in references 3 to 8, 11 
and 12. 

4.8 Computation of F Term in Inverse Problem 

It is pointed out in Chapter 3 that in the inverse or design problem the following condition of 
integrability should be used in the solution in the fluid flow along the S 2 stream filament. 

F • V x F = 0 (4.86) 


When general curvilinear coordinates are employed, substituting equations (4.12), (4.17) and (4.63) into 
equation (4.86) gives 


d 


d 

V 






(4.87) 


Then 


U 

f ) 
‘2 

+ 

rx 1 d 

v 

= 


l(xl) o S 




0 



(4.88) 


or 


(4.89) 


F 2 

F 2 ' 

+ If,, - 

\ 

F l/ \/ a ll 


i^vl 

J * , >" dx ! 
o 

V ) 


The procedure of calculating f 2 is similar to f r of references 2, 12 and 16. Compute f x from 
equation (4.68a), from (4.68c), the derivative of tjt* with respect to x 2 , then integrate f 2 along 
x 1 from (x 1 )^ where f 2 is known. It is convenient to choose (x 1 ) c at a station where Fr is equal to 

zero because F 2 can then be computed from F r 

4.9 Computation of Components of Normal and In Direct Problem 

In the direct or analysis problem, the known S 2 surface offers two relations between the n- or F 

components. Since n is normal to W, 


W = n,w + 


Nj W 1 ^2 W 2 , ^ ^ _ 0 

n,w 2 4 - n^Wp = -^=. + -i=. ~T=- + IN »s vv ^ 


Ja" fin 


/a 22 V a 22 


(4.90) 


or 


N l W* sin $ l2 + N 2 W 2 sin 0 12 + N^W^, - 0 


(4.91) 


then 


90 


I N 1 , N 2 2 I 

W„ = - _ W 1 + _ W 2 sin 0 12 

* N N 

V<P P 


This equation is similar to the equation appearing in section 3.3. The ratio of the N-components may be 
obtained by the following consideration: Along the intersection of the S 2 surface and a constant x 
surface 


J 11 fd^l / / * 

N v> Ux 1 ; Ux 1 J x 2 =(x 1) c 


Along the intersection of the S 2 surface and a constant x 1 surface 


N 2 _ _ r d<p ^22 __ _ r dip 
N,„ P, 2 p v 2 


dx V=(x 1 ) 0 


(y»22 sil 


Knowing these N-ratios can be computed from equation (4.92) and computed from W^. 

is then calculated by equation (4.68c) and F 2 by 


F ? = _ F, 
2 N * 
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CHAPTER 5 


TWO- AND THREE-DIMENSIONAL SUBSONIC FLOW IN TURBOMACHINES 

5.1 Series Expansion on Sj Surface of Revolution 

The intersection of an Sj surface with a mid-channel S 2 surface is a streamline lying in the mid- 
passage on the S t surface (fig. 5.1) and is referred to as the “mean streamline” in references 54 and 4 to 
6. In order to get an idea about fluid flow along the S 2>m filament with reference to the flow in the 
three-dimensional channel, the fluid flow along the mean streamline was examined in that reference based 
on the solution of compressible (subsonic) flow in a typical two-dimensional turbine cascade obtained by 
an electro-mechanical hand computer in the late forties. 

Figure 5.2 shows that as the inlet Mach number increases, the increase in the velocity along the mean 
streamline is larger toward the suction surface. The dip in velocity variation along the suction surface is 
also seen along the mean streamline. In figure 5.3 the shape of the mean streamline is compared to those 
of the mid-channel line and the mean chamber line. 

Figure 5.4 compares the variation of the specific mass flow along the mean streamline to the one- 
dimensional variation, the passage width. Inside the blade passage the former is a few percent greater 
than the latter, however the difference at the leading edge is large (25 percent). The effect of the blade 
thickness extends 40 percent upstream and downstream of the blade is also apparent. 
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The variation of tangentail velocity and its rate along the mean streamline are shown in Figure 5.5. The 

two begin to vary at 60 percent chord length upstream of the blade and end at 60 percent chord length 

* 

downstream of the blade. 

Based on the preceding observation a very simple, quick approximate solution of the inverse or direct 
problem for fluid flow along an arbitrary Sj stream filament or revolution was proposed in references 5 
and 6. Any flow variable on the mean streamline is considered a function of the meridional coordinate 

t. 


q = q[*.(p(0] 


(5.1) 


Its total variation with respect to ? is 


dq _ dq _j_ dq d<f> _ dq ^ tan fi dq 

dt ct d<{> dt ct r d<f> 


(5.2) 


From the continuity and dynamic equation the following equations for the first derivatives, with respect 
to the angular coordinate p of the velocity components, are 


aw, 

dp 


r 1 + 

d« 


tan /3 d(rr/»W, ) 
rp df 


+ (W, + 2wr) 



(5.3) 


*A variation of tangential velocity, similar to that shown in Fig. 5.5 but simpler, was used in the 
calculation of fluid flow on S 2 m surface in Refs. 10 and 11. 
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(5.6) 


The variation of a velocity component Wj or other flow variables in the circumferential direction 
is given by the Taylor series 


fl(<A) = Am + “ ^m) 


1 

aq 

A 


(d<j>\ 

m 

2 (d<f> 2 ) 


+ 


(5.7) 


The mass flow between the mean streamline and suction and pressure surfaces are, respectively, 
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Mg = tt pW|d^ 


(5.8) 


and 


Mp = rr f' Pp ^>Wjd^ 

^ m 


(5.9) 


The mass flow between the streamline selected in the mid-passage of a radial-flow impeller and the 
suction surface, is 65 percent of total mass flow (ref. 6). The streamline distribution and the variation of 
the meridional velocity along the mean streamline are shown in figures 5.7 and 5.8, respectively. It is 
interesting to see that the final value of obtained in the series (three terms) calculation is quite 
different from the constant value assumed in the design and compares well with the value obtained in the 
accurate direct matrix solution (ref. 8), 

This method has been improved and further developed in references 55 to 57. Recently a computer 
code employing fourth order Taylor series has been developed (ref. 58). In the case of a transonic turbine 
cascade the result obtained by using the fourth order series is generally better than the result obtained by 
using the second order series (fig. 5.9). It is also seen from the figure that the result agrees with the 
experimental data better than that obtained by the time-marching method reported in reference 35, 
especially in the low supersonic region. Yet the code can be put on a TI-59 pocket calculator. 

(Computing time on this calculator is 3 min per station). 
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5.2 Series Expanison In jp-Coordinate 

The circumferential extension of flow properties on a mean streamline by series expansion along the 
S x surface of revolution, described in the previous section, can readily extend to a circumferential 
extension of flow properties on a mean S 2 stream surface by series expansion in the three-dimensional 
blade channel. Because S x flow surface of revolution is not assumed here, the three-dimensional series 
expansion method yields a full-three-dimensional flow field, whereas the use of series expansion on a 
number of S x stream surfaces of revolution yields a quasi-three-dimensional flow field. 


Partial Derivatives of Flow Quantities in Nonorthogonal Curvilinear Coordinates: Flow properties at a 


point in the three-dimensional flow channel are calculated by Taylor series expansion similar to 
equation (5.7). The First order partial derivative of W 1 and p, with respect to ip can be obtained from 
the three-dimensional vorticity components, equation (4.17), for isentropic irrotational flow, is as follows 
(refs. 59, 61): 


J&22 cos ^12 ^ a 22 r 

’ dx 2 


aw 1 
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C = 
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dx l 


pW'fii 


sin 


'12 




1 



and 


dp _ p 

d<p (7 - l)h 


— - JL { (pW 1 ) 2 + (pW 2 f - 2(p W 1 ) (p W 2 ) cos e i2 + (pWp) 2 J 


dp dp 


The higher order derivatives can also be obtained from the basic equations. 


5.3 Forming Successive S 2 Surfaces by Progressing Circumferentially from S 2 m 


Similar to that in the two-dimensional series expansion method, the space coordinate p of the S 2 
surfaces can be obtained by the following series (ref. 59) 


A <p = <p-<p m = {*2 - 


r > 

dp 

+ - (* 2 - 


\ ** > 

2 

m 

«*; 


+ ... 


( 5 . 11 ) 


where p is the circumferential coordinate of any S 2 surface and ¥ 2 is the value of stream function of 
that S 2 surface. 
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In reference to figure 5.10 the integrating factor B of S 2 surface in section 4.5 is related to its 


inlet value by 


JL = AiP (5.12) 

Bj (A*>)i 


The inlet plane Zj is taken sufficiently far upstream of the blade row where the flow is uniform. 

The difference of stream function ¥ between two adjacent S 2 surfaces is the mass flow passing through 
the shaded area in figure 5.10, i.e., 


A*2 





(5.13) 


With the use of equations (5.12) and (5.13) the <p coordinate of S 2 surface can be determined step-by- 
step as shown in figure 5.11. The velocity component W* and p at the S , S ,... are then calculated 
by Taylor series expansion. 

The method described in the preceding paragraphs was programmed in Fortran IV and employed to 
investigate three-dimensional flow in a compressor stator and a turbine rotor. Less than 3 min were 
needed on a M150 computer (0.4 MIP). The turbine rotor is designed for a pressure ratio of 1.4, inlet 
absolute Mach number of 0.4, and rotor hub tip ratio of 0.66. 
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Figure 5.12 shows the projection of the streamlines on S 2>m surface on the meridional plane. The 
blade shape and the velocity distribution around the blade is shown in figure 5.13. The intersections of 

the Sj and S 2 surfaces, with the through-flow cross-sections, are shown in figure 5.14. 

The twist of the Si surface in the radial direction (fig. 5.14) is seen, generally, radially outward 
toward the suction surface. As this turbine rotor is designed for radially nonuniform work output, the 
radial twist of S, surface is relatively large, but it is still under 3 percent. The differences in the normal 
distance between the adjacent S, stream surface, both in the flow direction and in the circumferential 

direction of the twisted S x surfaces, and the untwisted Si surfaces of revolution will explain the 

difference in flow fields between the two cases. 

5.4 Coordinate Transformation and Direct Expansion Method 

For a more accurate determination of the coordinates of the S 2 surfaces and flow properties on 
these surfaces, the following coordinate transformation is made: 


{ = x l (r,z) 
n = x 2 (r,z) 

= iM r >M 


(5.16) 


The surfaces of revolution formed by revolving the x 1 - and x 2 -coordinate lines around the z-axis and the 
S 2 surfaces are the new coordinate surfaces (fig. 5.15). From the continuity and dynamic equations, the 
following first order partial derivatives are obtained (refs. 60 and 61): 
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where 
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In addition second order partial derivatives are also obtained. Then 
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(5.20) 


av- 


This approximate method of solving three-dimensional flow has been programmed in Fortran IV. The 
complete solution was computed on a UNIVAC-1100 computer in less than 1 min. 


The method was used to calculate the three-dimensional flow field in the Chinese Academy of 
Sciences’s (CAS) reseach compressor (ref. 62). The mass flow is 61 kg/s, the rotor tip M is 1.4, the 
stage total-pressure ratio is 1.5, the hub-tip ratio at stator inlet is 0.49, and the number of stator blades 
is 37. Projections of S 2>in and its streamline projection on the meridional plane are shown in Figure 5.16. 

The blade shape and distribution of velocity on the k = 7 coordinate surface obtained by this method 
and the stream surface extension method are compared in figure 5.17. The difference between the two is 
small. A comparison of Mach number distribution with that obtained by the S 1 /S 2 iterative solution 
(see section 5.9) are shown in figure 5.18. Except near the leading edge, the result obtained from the 
present method is close to that of the three-dimensional solution. The relative twist of S x surfaces at 
the suction surface (J-4) is shown in figure 5.19. Rm is the value of R on the S 2 m surface of the same 
Sj surface. (R is the local value to the point j = 4, k = 1, . , . 11.) It is seen that the largest AR/R 
occurs a short distance from the hub wall and the maximum difference between the present solution and 
the full- three-dimensional solution also occurs there. It may be noticed that the maximum relative twist 
is rather small, being less than 1 percent. 
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5.5 Direct Matrix Solution of Subsonic Flow Along S x Filament of Revolution 


The method of employing general nonorthogonal curvilinear coordinates and corresponding 
nonorthogonal velocity components to express the basic equations governing the fluid flow along Sj and 
S„ stream filaments were developed during the sixties (refs. 28 and 29). Since 1969, computer codes 
using general curvilinear coordinates have become available for flow along the S 2 (ref. 40) and Sj 
(ref. 64) stream filaments. 


When nonorthogonal curvilinear coordinates (x 1 , x 2 ) are selected on the flow surface of revolution 
(fig. 4.1) the principal equation (4.48) is to be solved. The metric tensor a Q ^ involved in the equation 
is calculated from the given geometry of the blading in cylindrical coordinate system (I, <p) by the 
following equations: 


a 22 — 


ii ) 2 + r 2 [^T 
lax 1 ) lax l J 

iiT + r 2 (^T 

k a x 2 J Ux 2 J 


dt dt . .2 d<p d<p 
a l2 = + r 


dx 1 dx 2 


dx 1 dx 2 


cos — 


*12 




’22 


(5.21) 


Discretization with the use of second order central numerical differentiation formula leads to a nine-point 
star in the two coordinates (fig. 4.7). The resulting N algebraic equations for the N interior grid 
points is written as 
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[M]{*} = {P} 


(5.22) 


The coefficient matrix [M] is a tridiagonal matrix. In order to take advantage of the aero elements in 
[M], one row of fictitious grid points was added upstream and downstream of the blade and the only 
coefficient along the diagonals was numbered and stored (ref. 63), thus reducing the storage requirement 
for the coefficient matrices. For a problem of m = 47, n = 11, (517 grid points), only 28 000 coefficients 
are stored* Therefore, 32 K internal storage is sufficient. 

As in references 10, 11, and 63, the coefficient matrix is factored to a lower triangular matrix [L] 
and an upper triangular matrix [u] which has elements along the diagonal equal to unity. Thus, 
equation (5.22) becomes 

(Ll |u] {*} = <P> (5 ' 2! 


{♦} is then obtained by a simple forward and backward substituting process as follows: 


Calculate {Q} from [L] {Q} = {P} 


(5.24) 


Caluclate {¥} from [u] {¥} — {Q} 


(5.25) 
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The preceding direct factorization technique was available in the forties. It was First utilized to solve 
the S flow in a gas turbine (ref. 10) Only 385 interior points were involved in that problem, and the 
method woe tested on an IBM CPCE and a UNIVAC. Nine digits were used on the former and 11 digits 
were used on the latter. The results agree up to the fifth digit. The residual left at any grid point b less 

than one in the eighth digit. 

As in references 10, 11, and 63 density is obtained by storing the E - $ table (eqs. 4.39 to 4.41) m 
the computer. During the iteration process for * and the density terms on the right side of the 
principal equation (4.47), the relaxation factor for density equal to unity is good for quick convergence, 
whereas, the relaxation factor for ¥ should be reduced from 1 to 0.5 when the inlet Mach number 
increases from low subsonic to high subsonic. Depending on the inlet Mach number, usually 10 to 20 
iterations are required to reach the following convergence criterion (ref. 63): 


P — < 10 -5 


* V+1 - * v < 10- 4 


r 


(5.26) 


This computer code has been widely used by a number of people on a number of different computers. 
On a Sun 4/110 (32-bit, 7-MIPS) computer, 20 iterations for an S x problem of 500 grid points were 

computed in less than 30 sec. 

This code has been modified slightly by adding a nearly orthogonal C-mesh sub-system around the 
leading edge of the blade (ref. 65) (fig. 5.20). The ^-values at the boundary of this sub-grid-system are 
taken from the result of the solution formerly obtained with the H-shape grid system. This additional 
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calculation improves the flow field around the leading edge (fig. 5.20). Not only the peak velocity 
(fig. 5.21), but also the stagnation point at the leading edge are correctly determined as indicated in 
figure (5.22). For the blade shown in figure 5.22, the positions of the stagnation point corresponding to a 
number of inlet angles are shown in figure 5.23. 

Alternately an H-C mesh, such as that shown in figure 5.24, may be used in the computation of Sj 
flow to obtain flow computation around the leading edge accurately. The H-C grid is obtained by a 
numerical solution of Poisson’s equation incorporating finer grids in areas of the physical plane where 
large flow gradient exists. 

The computer code of reference 43 has been slightly modified to obtain the flow field through a 
turbomachine blade row having splitter vanes or tandem blades (ref. 66). Division of mass flow by the 
splitter vane (fig. 5.25) and the outlet flow angle are determined by applying Kutta-Joukosky condition 
to the trailing edges of the main blade and the splitter vane. Figures 5.25 and 5.26 show some of the 
results obtained in an investigation of the effects of replacing one half of the main blade by splitter vane. 
Figure 5.26 reveals that by substituting seven splitter vanes, which have half of the chord length and the 
same shape as the main blade (except around the leading edge), the peak velocity at the leading edge of 
the main blade is reduced to about the same value in a cascade of half the solidity. The Mach number at 
the second peak is also reduced. On the other hand, the outlet angle is reduced only 1.2° compared to 
2.3° in cascade of half the solidity. 

Iterative SOR and Direct Matrix Solution of Subsonic Flow Along S 2 Filament 

The computer code for the solving of the principal equation (4.84) for fluid flow along an S 2 stream 
filament, employing nonorthogonal curvilinear coordinates and corresponding nonorthogonal velocity 
components, was first programmed by the use of iterative relaxation method (ref. 64). In order to get an 
accurate solution, at least six stations are placed inside the blade passage (from leading edge to trailing 
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edge of the blade). Because of the large number of grid points involved in multistage machines (fig. 2.4) 
a large internal storage capability is required. However, such storage was not available at that time, the 
line over-relaxation iterative method was used. Usually 20 to 30 iterations are required to reach the same 
convergence criterion set for Sj calculation. (In this way, 1200 grid points can be used on a 32K 
internal storage.) Recently the same principal of equation has been programmed using the direct matrix 

method of solution. 

Similar to the Si flow passing through the blade row having splitter vanes, this code has also been 
slightly modified for calculating S 2 flow for fan blades, downstream of which the flow is divided into the 
core compressor and the outer annulus. 

Examples of solutions obtained by this computer code will be given in the next section for three- 
dimensional flow calculation and also in chapters 6 and 7 for transonic flow. 

The S 3 computer code is programmed for the inverse problem, i.e., the variation of V^r along the 
S 2 surface is specified by the designer and the condition of integrability is incorporated in the code (see 
section 3.3). Equations (4.92) to (4.94) are also incorporated in the program so that it may be used to 
calculate V^r from W . Thus the same code can also be used for direct problem. 

5.6 Three-Dimensional Subsonic Flow in Turbomachine 


Three-Dimensional Flow and Quasi-Three-Dimensional Flow in Turbomachines 

In general the solution for three-dimensional flow in turbomachines involves the use of a number of 
S, and S, stream filaments. Flow through a three-dimensional flow passage formed by two adjacent 
blades can be considered to consist of a Urge nnmber of thin S, stream filaments (fig. 5.27) and the floe, 



in each filament can be taken to be the same as that on the Sj stream surface at mid-height of the 
filament (fig. 5.28). Similarly, the flow through the passage can be considered to be consist of a large 
number of S 2 stream filaments (see fig. 5.29) and the flow in each filament can be taken to be the same 
as that on the S 2 stream surface in the mid-circumferential position of the filament. The two families 
are related in the following manner: The stream surface of one family is formed by joining the 
corresponding streamlines lying on the stream surfaces of the other family. Hence iteration between the 
solutions of the flow through the two families of stream Filaments is necessary to obtain the three- 
dimensional flow. Theoretically the solution becomes exact when the number of filaments becomes very 
large. However, for engineering accuracy, 8 to 10 filaments in each family is sufficient. 

For axial-flow turbomachines the twist of the Sj stream surface is usually not large. In the so- 
called “quasi-three-dimensionaP solution all of the Sj-stream surfaces are assumed to be surfaces of 
revolution. In that case the radius coordinate and the normal thickness of the filament do not vary in the 
circumferential direction. Any variation in flow direction can be determined by the use of only one S 2 
surface. In the following, the S 2 surface lying in the mid-channel, the S 2>m surface, will be used 
(fig. 5.30). This choice is, of course, arbitrary. But, it is preferred, since (1) the shape of that surface 
bears some resemblance to the mean surface, (2) the flow on that surface is an approximate average flow 
in the channel, (3) the calculation of the flow on that surface is the first step of practical three- 
dimensional blade design process, and (4) calculation of flow on more S 2 surfaces will be involved in the 
full- three-dimensional solution. The geometry of the Sj stream surface determined by an axially 
symmetrical solution or by an averaging process of the flow in the circumferential direction would give a 
quasi-three-dimensional flow solution different from the one described herein. 
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Three-Dimensional Blade Design Procedure 

It is suggested that the following steps may be followed in a three-dimensional blade design process: 

(1) Determine flow path in the meridional plane. Simplified radial equilibrium calculation is made 
at a number of stations between rows (ref. 1), If desired, this flow path can be modified in the 
S 2 m calculation in the third step. 

(2) Estimate radial variations of stator stagnation pressure loss and rotor isentropic efficiency and 
compute the corresponding entropy increases across the blade row (see section 2.5). Assume linear 
variation of entropy in the flow direction except at the passage shock. It is desirable to check the 
loss with boundary layer calculation in later cycles. 

(3) Assume the variations of rjr- x and (V^r)/(V^r). along the mid-channel S 2 filament (see section 
5.1). Solve the S 2m principal equation for $ 2 ( section 4.6). Compute density and velocity 
components. Compute coordinates ^?(x*,x^) or (f,£>) and examine the shape of S 2 for its 
smoothness. 

(4) Obtain the geometry of the Sy stream filaments of revolution, [r(z),r(z)] from the meridional 
projection of streamlines on S 2 m . Design the blade section by the mean-streamline/series 
expansion method (section 5.1) or select standard blade section on the conical surface that 
approximates the S x surface in the blade passage. Solve the S x principal equation for ¥ 

(section 4.4). Compute velocity components and density. If the p or W distribution around 
the blade has to be improved, repeat the design or blade selection process in step 4. An inverse 
calculation is available to modify the blade shape for a better W or p distribution around the 
blade (ref. 67). The new values of V^r and r obtained are inputs to the next S 2 m calculation. 
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(5) Repeat step 3 with the (V^r) and r values obtained in step 4. Modify the passage geometry in 
the meridional plane if desirable. Modify the stacking line in the range permitted by stress 
consideration if desirable. 


(6) Steps 4 and 5 are repeated until acceptable convergence and desirable velocity and pressure 
distributions around the blade surface and along the hub and casing walls are obtained (see 
fig. 5.31). 


The S| surface near the hub or casing wall is usually assumed to be surface of revolution. If the 
angle between the Sj surface there and the hub or casing wall in the meridional plane is denoted by a 
(fig. 5.6), and if the Sj surface there is a surface of revolution, then 


d<p 


= 0 


By using the following relation between a and the velocity components 


(5.27) 


tan a = 


Wr 

Wz 


and the following equation of vorticity component, 


(V 


i a w . 

x v), = I 1j 

r dp 


1*V) 

r dz 


(5.28) 


(5.29) 


there is obtained (ref. 68) 
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1 do_ = 1 1 d(V) _ p v yj (5.30) 

t r d<p W| r dn 

Equation (5.30) shows how (V^r) and the vorticity affects the twist of the Sj surface. In particular if 
Sj surface is a surface of revolution, da/d<p = 0 and equation (5.30) becomes 


d(Vgr) 

dn 


= r(V X V), 


(5.31) 


If the (-vorticity component is equal to zero, then 


d( V 0 r ) _ Q (5.32) 

dn 

Thus, in the third step of the design process described in the preceding section, satisfaction of 
equation (5.32) in the values of (V^r) specified near the hub and casing wall would help keep the S 1 
stream surface near the wall being surface of revolution. 

Determination of the Circumferentially Uniform Conditions Par Upstream and Far Downstream of the 
Blade 

In the solution of S x flow the uniform condition (in the circumferential direction) at stations far 
upstream and far downstream of the blade row is required. This is true not only for all of the blade rows 
in a multistage machine (fig. 2.4), but also for the rotor and stator in a single stage machine (fig. 5.32). 
The far field boundary condition upstream of blade row S in figure 5.31 is obtained as follows: 

(1) Calculate S 2 m flow with all of the blade rows in the flow passage. 
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(2) Take off blade row R and make another S 2m flow calculation from station AA to station BB. 
In this calculation the $ distribution obtained at station AA in step 1 is used as the boundary 
value at that station. At the far upstream station BB, a boundary condition of the second type 
(i.e., Vr = 0) is considered. 

(3) With the ¥ distribution at BB and a second type boundary condition (i.e., Vr = 0) specified at 
far downstream station CC, an S 2 m calculation from station BB to station CC is made. 

(4) The fluid state obtained at station BB in step 2, and the streamline shape and the fluid state 
between station BB and CC obtained in step 3 provides the geometry of the Sj surface and the 
boundary conditions at the two stations for the calculation of flow along a number of S L stream 
filaments. 

The far downstream condition for blade row R is obtained in a similar manner as indicated in the 
figure. 


5.7 Three-Dimensional Flow in a High Subsonic Compressor Stator 
Designed by Quasi-Three-Dimensional Flow Method 

The quasi-three-dimensional blade design method described in the preceding paragraphs has been 
used to design a high subsonic stator of a single stage transonic fan. The projection of the stator blade 
on the meridional plane is shown in figure 5.33. The inlet Mach number at the hub is 0.8. After 
completion of the blade design a three-dimensional analysis was immediately carried out. The quasi- 
three-dimensional design took five cycles of iteration, whereas the full-three-dimensional analysis took 
only three more cycles of iteration. 
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It was customary in the past that entropy gradient was considered in the S 2 calculation (fig. 2.5), 
whereas isentropic assumption was considered in subsonic Sj calculation. Iteration between the two 
immediately shows divergence of the solution (fig. 5.34). When the same value of entropy was used at 
the same point on the two surfaces, convergence was quickly reached after only two more cycles of 
calculation. The criterion used for both convergence consistency is 


M_(n) ~ M o(") 


M o(») 


< e 


(5.33) 


e was reduced to less than 2 percent. Figure 5.35 shows the streamlines obtained at the root in Sj and 
S, m calculation in the last cycle. 

it jUl 

In the followup, full-three-dimensional analysis calculation the Sj surface is no longer considered as 
surface of revolution. The nonorthogonal curvilinear coordinates are still placed on the surface. The 
principal equation is the same as equation (4.84) except for the last term in M. Instead of 
(2wr p sin a sin 2 0 l2 ) it is now (2 <A/> sin 2 0 12 ). Actually the computer code needs only one S t code and 
the general S x code. By putting w 3 = u sin a, the same code becomes the surface of revolution 

code. 


Due to the twist in the general S x stream surface, it is not possible to impose the periodic condition 
immediately downstream of the trailing edge of the blade. Instead the boundary condition of the first 
kind is used there, that is, ¥ is specified there. After a certain distance downstream the periodic 
condition is again imposed. 
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The difference between the two solutions is indicated by the twist of the Sj surface, that is the 
deviation of Sj surface from surface of revolution. Figure 5.36 shows the amount of deviation of the Sj 
surface number 3 and number 5 of figure 5.37 at a number of stations along the flow direction. It is 
apparent that the deviation of surface number 3, located nearly at mid-span, is larger than that of 
number 5. The deviation measured by the difference in radius, Ar, is less than 1 percent of the radius r. 
Figure 5.37 also shows twist of S x stream surfaces near the leading edge and the trailing edge of the 
stator. 

The difference in the streamline distribution, the variation of V^r along the intersecting streamline 
of Sj and S 2 surfaces and the variation of the angular thickness of S 2 m filament as determined in the 
Sj calculation, are shown in figures 5.38 to 5.40, respectively, and are small. Yet these relatively small 
differences in radius or filament thickness causes a relatively large difference on the Mach number 
distribution over the blade surface (fig. 5.41). This phenomena can be explained by examining the 
difference in the variation of r in the two cases. Figure 5.37 shows that r near the suction surface at 
the blade leading edge is smaller in the full-three-dimensional flow than in the quasi-three-dimensional 
flow. Thus, the peak Mach number on the suction surface obtained in the full-three-dimensional analysis 
solution is higher than that obtained in the quasi-three-dimensional design calculation. It exceeds the 
sonic velocity slightly. 

5.8 Three-Dimensional Flow in CAS Research Compressor (Subsonic Case) 

CAS Research Compressor 

The research compressor designed for experimental investigation of three-dimensional subsonic and 
transonic flow in an axial-flow compressor is a single-stage compressor with inlet guide vanes. It is 
designed for a stagnation pressure ratio of 1.5 and isentropic efficiency of 0.85. The tip speed is 400 m/s, 
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the inlet hub-tip ratio is 0.4. The corrected mass flow per unit frontal area is 188 kg/m 2 s and the 
average inlet Mach number is 0.616. 


Figure 5.42 shows the flow path with the intersecting lines between the meridional plane and Sj 
surfaces. The hub and casing contours are so chosen as to reduce, as much as possible, the sever flow 
conditions at the tip of the rotor and at the hub of the rotor and stator. The relative Mach number of 
air entering the rotor varies from 0.68 at the hub to 1.34 at the tip. The MCA blade shape (laid out on 
conical surfaces) was used for the tip elements. The DCA blade shape was used for the elements in the 
mid-span and hub regions. Blade solidity varies from 2.3 at the hub to 1.3 at the tip, and the diffusion 
factor at the tip is limited to less than 0.4. The absolute Mach number of air entering the stator varies 
from 0.80 at the hub to 0.58 at the casing. Solidity of the stator blade varies from 2.0 at the hub to 1.1 
at the casing and the diffusion factor is limited to 0.5 at the hub. The turning angle at the hub is 
slightly less than 45 degrees for the rotor and slightly less than 45 degrees for the stator. Aerodynamic 
design details are given in reference 62. 

Overall Performance 

Rotor test data were taken over a range of mass flows from maximum flow to near stall condition at 
60, 70, 80, 90, and 100 percent of equivalent design speed. For each mass flow, measurements were taken 
at seven radial positions. Measured outlet stagnation pressures, stagnation temperatures, stream static 
pressure and flow angles were corrected for Mach number effect. 


Overall stagnation pressure and stagnation temperature ratios were obtained from a mass average of 
the survey data at the rotor outlet and the pressure and temperature measured at the inlet station. 
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Strain-gauge type and Setra capacitance type transducers were used in measuring pressures. Copper 
resistance temperature detectors were used in measuring temperature at the inlet station and nickel 
chromel silicon-nickel thermocouplers were used in measuring temperature at the outlet station. Flow 
through the rotor was determined by the calibrated bell-mouth inlet. 

Rotor speed was measured by a magnetic pick-up in conjunction with a gear mounted on the drive 
motor shaft. A phase difference type torque-meter was used to measure the torque of the drive shaft. A 
3000 kW D.C. motor and gear boxes were used to obtain speeds up to 12 000 rpm for the research 
compressor. 

The estimated accuracy of these conventional measurements is as follows: 


Inlet pressure, mm H 2 0 ±10 

Outlet pressure, mm H 2 0 ±25 

Temperature, K ±1 

Mass flow, % ±0.5 

Speed, % ±0.03 

Flow angle, deg. ±1.5 

Torque, (0-400 kg/m), % ±1 


For overall performance L2F measurements were only made at the inlet and outlet station. At the 
inlet station, data were taken only at one circumferential position, whereas at the outlet station, data 
were taken at 16 circumferential positions. 
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The accuracy 


of the L2F measurements is estimated to be 


Flow speed, % 
Flow angle, deg. 


±1 

±1 


The overall performance of the rotor is plotted in figure 5.43. The plotted data present rotor 
stagnation pressure ratio and isentropic efficiency as a function of equivalent mass flow for rotative 
speeds of 60, 70, 80, 90, and 100 percent of design speed. 

Detailed L2F Measurements 

Detailed L2F measurements were first made for the peak efficiency point of 70 percent of design 
speed. The equivalent mass flow was 144.8 kg/m 2 s, the stagnation pressure ratio was 1.21, and the 
isentropic efficiency was 93.2 percent. Flow speed and flow angle were recorded at 8 axial stations, 6 to 
14 radial positions, and 10 to 16 circumferential positions. All of the measured points are on the S 2 m 
and Sj stream surfaces obtained in the design calculation. 


Three-Dimensional Flow Calculation 

Before the result of the three-dimensional flow field calculation can be properly compared with the 
L2F measurements, suitable flow parameters measured at the outlet station should be used in the 
calculation. The probe measured outlet stagnation pressure p 2 ° and laser measured outlet flow angle 
« 2 , shown in figure 5.44 as functions of radius are chosen. The main input to the three-dimensional 
analysis calculation is then: flow path (in meridional plan), blade geometry (coordinates given on design 
conical surface or on manufacturing templates), inlet pressure and temperature, rotor speed, mass flow, 
flow coefficients, and outlet stagnation pressure and outlet absolute flow angle. 
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The first step in the three-dimensional calculation is computation of flow along an S 2 jn filament 
about mid- way between two adjacent blades (fig. 5.31). The radial distribution of rotor efficiency 
determined by temperature measurement, is shown in figure 5.45. The measured value of torque is 
considered to be more accurate than the measured value of a small increase in temperature across the 
rotor. The temperature determined radial distribution of rotor efficiency, used in the initial calculation, 
has to be corrected so that the radial distribution of stagnation pressure in the calculation at the outlet 
station should be at least approximately equal to that obtained by the pitot tube measurement. The 
span-wise averaged value of the corrected efficiencies has to agree with the rotor efficiency determined by 
torque measurement. The resulting distribution of rotor efficiency also determines the radial distribution 
of V^r at the outlet, which is compared to the values used in the design (fig. 5.45), i.e., 100-percent 
design speed. The comparison shows that at the peak efficiency point of 70-percent design speed, the load 
distribution along blade span is quite different from that at the design point. The load near the tip and 
the efficiency distribution shows that the efficiency at the hub is lower than that at the design point. It 
seems that is caused by the appearance of a higher angle of attack in the hub region. 

In the initial calculation the variation of V^r along the flow direction is taken to be the same as 
that used in reference 11. The variaiton of the stream filament thickness is estimated according to the 
blade thickness distribution. (In a later calcualtion, however, both estimates are not need, because their 
variations are given by S x solutions.) 

With these input values, solution of S 2 m flow is obtained. The axial velocity at the outlet station 
and the experimental flow angle a 2 yield a new set of values for V^r at the outlet station. The S 2 m 
calculation is repeated with this set of values until convergence of the solution is obtained. 
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The converged value of outlet V^r is then used with the experimental value of outlet stagnation 
pressure to correct the radial variation of rotor efficiency. The modified values are used to repeat the 
preceding calculation of S 2 m until the rotor efficiency value converges. 

From this solution of S 2 m the geometry of 11 surfaces of revolution are formed and interpolated 
values of coordinates of the corresponding blade surfaces and the variation of the S 1 stream filaments 
are obtained* 

After these solutions of the flow on 11 Sj surfaces are obtained, a new S 2 m surface, together with 
the corresponding variations of V^r and stream filament thickness, is obtained. The flow on this new 
S 2 m surface is computed in the same manner as previously discussed. 

The whole calculation process is very similar to the quasi-three-dimensional calculation described in 
section 5.7, except that input values to the whole calculation are different. The solution of the complete 
procedure converges rapidly (see fig. 5.46) and only four computation cycles are actually required in the 

present case. 

This quasi-three-dimensional computation could be conveniently followed by a full three-dimensional 
computation as suggested in the preceding section. It was found in the present case that the difference 
between the two solutions is rather small. The twist of the Sj stream surface obtained is very small. 
Comparisons were made only between the experimental measurement and the quasi-three-dimensional 

solution. 
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The following sketch shows the various steps invloved in the present calculation: 

Input Data 



Comparison and Discussion 

From the calculation procedure presented in the preceding section it is seen that the result of the 
theoretical calculation depends on the accuracy of the experimental data put into the calculation. 

Because the measured p 2 ° and are believed to be relatively more accurate than the others, they are 
chosen as the input for the calculation, with only a small adjustment in the values (see fig. 5.44). 


On the other hand, errors in temperature increase across a single rotor could be quite large, 
therefore, a relatively large adjustment in local rotor efficiency was made for the overall rotor efficiency 
determined by the torque-meter to be satisfied (see fig. 5.45). 
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Figure 5.47 shows the radial variation of the relative flow angle at the inlet and outlet station 
obtained in the theoretical calculation for 70-percent design speed. It is seen that their variations are 
quite similar to those of the design values. 

The corresponding radial variations of relative Mach number at the inlet and outlet station are 
shown in figure 5.48. Although the inlet Mach number is only slightly greater than 0.8 at the tip, the 
calculation shows that supersonic velocity occurs near the leading edge of the blade in the tip region. 


The absolute velocities of air flowing past the rotor blade channel obtained in the L2F measurement 
and in the theoretical calculation on the S 2 m surface are compared in figure 5.49. In general, the two 
agree closely in the trend of variation and agree reasonably well in magnitude. Toward the exit, the 
calculated value is lower than the L2F measured value. Perhaps this is caused by the inadequacy in the 
magnitude of the mass flow coefficient used in the theoretical calculation to account for the effect of an 
annulus boundary layer. 

The variation of relative flow angle obtained by the L2F measurement along a streamline near the 
blade suction surface, along the mean streamline and along a streamline near the blade pressure surface, 
all on an S t surface at 40-percent blade height, are compared to the calculated value in figure 5.50. In 
general, the agreement is pretty good. Similar to that noticed in the design calculation, the air has an 
overturning in the blade channel at the outer radius even at 70-percent design speed. Figure 5.51 shows 
constant Mach contours on three S ^ surfaces. 
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CHAPTER 6 


TRANSONIC FLOW ALONG Sj AND S 2 STREAM FILAMENTS 

Intensive research work on the development of a relatively simple and quick method of computing 
the transonic flow along Sj and S 2 relative stream surfaces has been carried out during the last 5 years 
at the Institute of Engineering Thermophysics, CAS. Once the quick methods for solving transonic Sj 
and S 2 flow have been completed, they can be readily utilized to obtain a quick solution of a three- 
dimensional transonic flow. 

Because the computer time required to solve a steady flow problem by using a time-dependent 
approach is quite long, a solution for the steady flow problem using much shorter computer time has been 
sought through a number of approaches. The mean streamline/series expansion method, as seen in 
sections 5.1 to 5.4, is quite simple and is applicable to a low transonic flow without strong passage shock. 
When strong passage shock in the the blade passage is considered, it immediately appears that using 
existing methods of calculating the supersonic flow and subsonic flow, separately, on the two sides of the 
passage shock should work. For instance the method of characteristics can be used for the supersonic 
region while the direct matrix solution can be used for the subsonic region. 

6.1 Transonic Flow Along S l Stream Filament of Revolution Solved by Separate 
Region Computation with Shock Fitting 

Calculation of Supersonic Inlet Region 

In figure 6.1 a typical compressor cascade with a supersonic inlet flow is shown with the bow- 
wave/passage-shock. The bow wave region upstream of bow wave can be easily calculated by the method 
of characteristics. References 37 and 38 give detailed descriptions of this method of solving supersonic 
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flow along an arbitrary stream filament of revolution. The nonuniform increase of entropy downstream of 
the detached bow wave is included in this computation. The calculation is carried out from one blade to 
the next in the circumferential direction for several blades until the periodicity condition is satisfied. 
Among others, the theoretical inlet flow angle which corresponds to the inlet Mach number Mj, is 
determined. It was noticed in the calculations that the periodicity of velocity can be satisfied easier than 
that of entropy (ref. 37). (Usually the former occurs at the fourth blade channel while the latter occurs 
at the seventh blade.) 

Location and Shape of Detached Bow Wave 

In the past, the location and shape of the bow wave were usually taken from empirical data of an 
isolated plane or axially symmetrical body (refs. 70 to 72). Based on experimental data given in 
reference 73, modification of such data for plane cascade has been done in reference 74. For three- 
dimensional applications the location and shape of the bow wave on the Sj surface of revolution is 
required and their determination by theoretical calculation has been presented in reference 75. 


In calculating the transonic flow field for a cascade test, the location and shape of the bow wave may 
be taken from the middle line of the bow shock band obtained through optical measurement. For the 
design of a transonic cascade, they may be taken from empirical data (refs. 70, 71, 73, and 74) or may be 
taken from the shock capturing solution described in section 6.2. 


Computation of Subsonic or Transonic Region Downstream of Passage Shock 

Depending on the inlet Mach number, the flow downstream of the passage shock may either be 
subsonic or transonic. The subsonic Sj code (ref. 63) is modified to take passage shock as the inlet 
boundary and to take the stream function variation immediately downstream from the passage shock as 
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the Dirichlet boundary value at the inlet station. The nonorthogonal curvilinear coordinate system is 
especially suited for use in the downstream region (see fig. 6.2). 


The streamline distribution and constant Mach number contours calculated for DC A 2-8-10 cascade 
(ref. 73) at an inlet Mach number of 1.11 and a static pressure ratio of 1.35 are shown in figure 6.3. 

In figure 6.4 the variations of three important flow variables V^r, 0 , and M, along the three 
streamlines are shown. Notice that the change in the flow direction, or angle /9, is large at the curved 
part of the passage shock, whereas the jumps in M and V ^r are large near the suction surface. 

Similar data obtained on a number of Sj surfaces gives the shape of the mean (and other) S 2 
stream surface, the thickness variation of the corresponding S 2 stream filaments, and the variation of 
V^r over the S 2 surface. The accuracy of the calculated flow along the S 2 surface depends on the 
accuracy of this data. 

Improvement of Location and Shape of Passage Shock 

According to the calculation for supersonic flow upstream of the passage shock and the assumed 
location and shape of the passage shock, the Mach number, just upstream of the shock, the angle included 
by the shock, and incoming flow are completely determined. Then the turning angle across the shock, the 
Mach number just downstream of the shock, the pressure increase across the shock, the density ratio and 
the entropy change, etc. are calculated by the Rankine-Hugoniot relations. If the assumed location and 
shape of the passage shock are correct, these calculated values of aerothermodynamic quantities just 
downstream of the shock and those at the inlet boundary obtained from the solution of the subsonic flow 
field, must be the same. Based on this criterion, the correct location and shape of the shock can be 
obtained. This adjustment process is incorporated in the computer code of reference 38. Figure 6.5 
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shows the agreement in W and /} in the case of a DCA 2-8-10 cascade calculation. On the UNIVAC 
1100, computation time for a supersonic region with 90-by-20 mesh points is less than 2 min, and for 
subsonic region with 13-by-40 mesh points, computation time for each iteration is less than 1 min. To 
obtain the convergent result, the total computation time is about 10 min. In reference 76, this method 
was applied to compute the transonic flow along three Sj stream filaments in the DFVLR rotor and 
compare the flow with experimental data. 

In case the flow downstream of the passage shock is transonic, one of the shock capturing code 
described in section 6.2, may be used. 


6.2 Transonic Stream-Function Principal Equation Solved With the 
Use of Artificial Compressibility 

In the treatment of a transonic flow calculation, the use of the stream function principal equation is 
preferred to the use of the potential function principal equation. Many good methods are available for 
quick, accurate solution of the elliptic principal equation corresponding to the direct problem of subsonic 
flow along an S x surface. For transonic flow, however, the stream function approach has inherent 
difficulty since gas density is not uniquely determined from the mass flux obtained from the derivatives of 
stream function $ after each cycle of iteration of the ¥ field. This is why long computer time was 
spent to obtain the steady solution through artificial unsteady calculations. Now, with the method of 
introducing artificial density through the density term (refs. 76 to 78, and 40), the stream function 
principal equation is elliptic throughout the whole flow field, and quick methods are available for solving 
the elliptic equation. Several attempts (refs. 41, 42, and 79) to assess this approach are described briefly 
in the following section. 
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Principal Equation 

The following conservative form of the stream function principal equation is obtained by substituting 
equation (4.46) into the dynamic equation (4.33): 


dx 1 


Ai- 


i a* A lat 

. - a 2 

P dx 2 


P dx x J 


( 6 . 1 ) 


d 

dx 1 


1 a* 1 a* 

a 2 A 3 

p ax 2 




where 



It may be noticed that in the stream function solution it is easy to include entropy variation in the 
flow field. The density calculation includes entropy increase as follows: 
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p = p 1 
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( 6 . 2 ) 
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The principal equation (6.1) also includes a variation of entropy with respect to x 2 . Entropy increase 
across a shock is calculated by the usual formula. 

Method of Solution 

The principal equation (6.1) is a mixed-type equation. It can be solved by modifying the density so 
numerical dissipation necessary is introduced in the supersonic region. For the present treatment, in 
which general nonorthogonal curvilinear coordinates are used, modified density is selected as follows: 

P = P ~ /*/> 8 As ~ p - t* 


W A W 2 A 2 


( 6 . 3 ) 


where 
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C = 0.5 - 3 


Determination of Density 

In transonic flow, the density is not uniquely determined by the value of mass flux, through the E-$ 
method usually used in subsonic flow cannot be used. Density is now determined in the following 
manner: By the use of equation (4.46) and artificial density, the dynamic equation (4.33) can be written 
as (see eq. 4.41) 
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(6.5) 


This velocity gradient equation is now used to obtain W 2 along the x 2 coordinate line by 
integrating from an initial x 1 coordinate line. W 1 and W 2 at grid points on this line are first 
calculated using equation (4.46), in which the density p is replaced by the artificial density p. Then 
W 1 on other x 1 coordinate lines are obtained by integrating from the initial line either in one direction 
(ref. 41), successively in reversed direction (ref. 42), or in positive and negative directions (ref. 79). 
Knowing W 1 , W 2 is calculated by 


wVw 1 = - 



(6.6a) 


W is calculated by 
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1 

w = [(w 1 ) 2 + (w 2 ) 2 + 2W*W 2 cos 0 12 ] 2 


(6.6b) 


and p is calculated by equation (4.2). 

Discretized algebraic equations obtained by the principal equation (6.1) by conservative central 
difference formula over a whole flow field can then be solved by a suitable method for the solution of an 
elliptic equation. In the direct matrix method of solution the coefficient matrix [M] of the discretized 
equation (6.1) is now a function of p. To prevent decomposition of the coefficient matrix [M] in every 
cycle of iteration, equation (5.0) is transformed into (ref. 41). 

[M'] [A V’] (n+1 ) = p [Rj (n ) ( 6 


so that the decomposition of [M'j is needed only the first time, and subsequent solutions can be obtained 
with relatively few forward and backward substitutions. 

In addition to this direct matrix solution, strong implicit approximate factorization procedures 
(refe. 80 and 81) and vertical line relaxation procedures were also programmed in reference 41. 

It was found when calculating for a transonic cascade with a mesh of 61-by-ll, the strong implicit 
algorithm and the relaxation algorithm took about 2 sec per iteration on an UNI VAC- 1 100 computer 
whereas the direct matrix algorithm took about 3 sec on the same computer. But convergence obtained 
after 40 iterations (120 sec) in the latter is better than the convergence obtained after 100 iterations 
(200 sec) in the former. The Mach number obtained by the direct matrix algorithm for a turbine cascade 
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and a compressor cascade are compared with other theoretical calculations and experimental data in 
figures 6.6 and 6.7 respectively. Fifty to sixty and 11 to 15 grid points were used, respectively, along x 1 
and x 2 coordinate lines, and 200 to 300 iterations were carried out to reach an ^ eve * 10 3 * ^ * 8 

interesting to see that in the case of the turbine cascade, the result obtained by the present method 
agrees quite well with experimental value. It seems that both this method and the mean 
streamline/series expansion method gives a Mach number distribution in the supersonic region better 
than that given by the time-marching method. In the case of the T-l compressor cascade, it is seen that 
the present method captured the second shock better than the relaxation method and time-marching 
method, and that the position of the captured first shock is a little ahead of that given by the experiment 
data. 

Improvement of Velocity Distribution for Transonic Blading 

On the basis of this treatment of the transonic direct problem, a method of solution for the transonic 
inverse problem was developed and reported in reference 67. This method is very useful in improving the 
aerodynamic performance of cascade blades. For instance, the velocity distribution on the suction of 
T x (18^1^)08 cascade shows two high velocity peaks (fig. 6.8). This method enables one to modify the 
velocity distribution so that the magnitude of the two velocity peaks is reduced, but the blade circulation 
is kept the same (see fig. 6.8). Blade coordinates obtained in the inverse solution are shown in figure 6.9. 
It is seen that the change in the blade shape is very small. Maximum blade thickness and outlet flow 
angle are practically the same. 

6.3 Effect of Axial-Velocity Density Ratio and Viscous Effect 

In comparing theoretical calculations presented in the preceding sections with experimental data, it 
is important to notice that the cascade experimental data is usually obtained at a certain axial-velocity 
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density ratio. In the case of compressor cascade it exceeds unity. Therefore, the theoretical calculation 
should consider a variation of the stream filament r corresponding to the axial velocity density ratio in 
the experiment. 


It is also important to consider the viscous effect. If flow does not separate from the blade, the 
viscous effect is reflected by the reduction in flow area due to boundary layer development on the blade 
surface and the increase of entropy in the fluid. 

At a relatively low supersonic inlet Mach number, inclusion of the effect of entropy increase in the 
theoretical calculation (shock capturing method) seems to give a pressure distribution in agreement with 
the experimental data (ref. 73) (fig. 6.10). The variations of V tf r, M, and fi along three streamlines 
are shown in Figure 6.11. 

In the case of an inlet Mach number of 1.11, the inclusion of both entropy increase and axial** velocity 
density ratio (fl = 1.16) in the theoretical calculation yields a fair comparison between the calculated and 
experimental pressure distribution (ref. 73) around the blade surface (fig. 6.12). The variations of V^r, 
M, and ft along three streamlines are also shown in figure 6.13. 

At an even higher inlet Mach number of 1.34, probably due to boundary layer separation, an fl 
value of 1.18, which is much higher than the experimental data of 1.04, is required to bring the pressure 
distribution close to the experimental level (fig. 6.14). The corresponding variation along three 


streamlines is shown in figure 6.15. 



Calculations (ref. 85) have also been made to obtain the flow along three S l surfaces at 68-, 89-, 
100-percent blade height of the DFVLR rotor. Here an adequate knowledge of the relatively large 
variation of r is very important for results to be close to the experimental data. 


6.4 Transonic Flow Along General Stream Filament 

In a general S l stream filament the radial coordinate of the central Sj surface and thickness of the 
filament vary in the flow direction and circumferential direction (see fig. 6.16). In transonic flow, there 
are abrupt changes in r and r in the shock region. The method for solving transonic flow along an S x 
stream filament of revolution given in reference 41, has been extended in reference 86 to that along an 
arbitrarily twisted S l stream filament. 


In a full-three-dimensional calculation, the shape of a general Sj stream surface is obtained by 
joining corresponding stream lines on a number of S 2 stream surfaces. In reference 87, general x 1 and 
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x 2 coordinates are placed on the central Sj surface of an S j stream filament. The third coordinate x 
is chosen in the radial direction (fig. 6.17). Then 
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where 


dr 
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( 6 . 9 ) 


The stream function principal equation in conservative form is 
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These equations are only slightly different from those on an surface of revolution. (Compare 
eq. (6.10) with eq. (6.1)). The procedure of solution for the transonic flow along a general S x stream 
filament is the same as that along an Sj stream filament of revolution (ref. 41). 

6.5 Transonic Flow on S 2 Stream Filament Solved by Separate Region Computation with Shock Fitting 

Similar to transonic S x surface flow, transonic S 2 surface flow can also be solved by a separate 
region calculation. However, the shock on an S 2 surface is different from that on an Sj surface. The 
shock usually begins at the outer casing and terminates a distance from the inner hub. Thus, the flow 
field is divided into three regions (ref. 39) (fig. 6.18). The first region AAjBjBA, is a supersonic flow 
region. It extends from the starting line AA 1? which may or may not coincide with the blade leading 
edge, to shock line BBj. Supersonic flow in this region is solved by the method of characteristics. The 
slopes of the two families of characteristics are 
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where i = 1,2 (when i = 1(2) the sign is + (-)). 


The consistency relation along the i family of characteristics is 
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The second region A 1 A 2 B 2 B 1 A 1 is a transonic flow region. The flow in this region is best 
calculated by a Taylor series, using the following partial derivatives with respect to x 2 : 
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The coefficients C lp C 2 ., C 3 ., E lp E 2 ., F 1;) F 2 . and G are given in reference 39. 

The third region BBjAjjAgCjCB is a subsonic flow region, the flow of which is solved by the SOR 


method. 


This method has been used in reference 39 to obtain the transonic S 2 flow in the CAS rotor and the 
rotor reported in reference 86. The computation was done on a UNIVAC 1100; CPU time was about 
3 min. For the high pressure (1.94) rotor of reference 86 meridional projection of the characteristic 
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network on S 2 and the mesh used for the subsonic region are shown in figure 6.19. Figure 6.20 shows 
meridional projection of the streamlines. The abrupt change of the streamline across the shock is clearly 
seen. The Mach number variation in the three flow regions is seen in figure 6.21, 

The variation of flow along the mid-channel S 2 surface can be seen more clearly by following the 
flow along three streamlines, namely, the one along the hub, the one at midspan and the one along the 
casing. Corresponding to the abrupt change across the shock (fig. 6.20) the decrease in Sj stream 
filament thickness is large (fig. 6.22) whereas the change in the relative flow angle (fig. 6.23) is small. 
Across the shock, the drop in relative Mach number (fig. 6.24) and absolute Mach number (fig. 6.26) and 
the increase in absolute flow angle (fig. 6.25), tangential velocity, (fig. 6.27) and pressure ratio (fig. 6.28) 
is relatively large. It is noticed from these figures that there is quite a large overturning of the air and 
overshooting of the tangential velocity below mid-span and near the trailing edge of the blade, 

6.6 Shock Embedding Elliptic Solution for Inverse Problem of Transonic S 2 Flow 

One of the findings in reference 3 has been generally accepted since it was published. That finding 
states that in the case of an inverse problem where the tangential velocity component is prescribed by the 
designer, the partial differential equation governing fluid flow remains elliptic as long as the meridional 
velocity component is lower than the speed of sound, even when the flow relative to the rotating blade is 
higher than the speed of sound. Indeed a great number of transonic turbomachines have been designed on 
that basis. However, if in the design of transonic machines, the design parameter selected, for instance 
V^r, is prescribed in the same manner as in the subsonic machines (fig. 5.5), flow discontinuity at the 
passage shock is not obtained in the solution. This kind of calculation is referred to as Level I 
calculation. If the values of V^r, s, and r prescribed in the design calculation are smooth, but have 
steep gradients at the passage shock, then the solution obtained indicates a passage shock. This is 
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referred to as Level II calculation. Level III calculation has the following characteristics: (a) prescribed 
values have proper discontinuity across the shock provided by transonic Sj calculations, (b) shock 
relations are embedded in the S 2 calculation, (c) S 2 solution obtained shows a clear shock discontinuity, 
and (d) one of the x 2 -coordinate lines coincides with the shock line (fig. 6.29). 

Embedding Shock Relations 

In the calculations of the direct problem of stream surfaces, usually the plane shock relations are 
used. This is equivalent to assuming that the shock surface is perpendicular to the stream surface. The 

error involved appears to be somewhat different for Sj and S 2 flows and its magnitude depends on 

particular configuration. In reference 88, this assumption is not made and fully three-dimensional 
relations are used across the shock. First, two stations, a very short distance apart, are placed along the 
shock line one immediately upstream and the other immediately downstream. The dynamic equation in 
the direction tangent to the S 2 surface and also to the shock surface is 

yj (+) _ yj (~) (6*17) 

ts ts 


The dynamic equation in the direction normal to the shock surface yields 
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The dynamic equation in the direction tangent to the shock, but not tangent to the S 2 surface yields 
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s(+) _ _ 7 ln [ 2 

““ 1 ~ = ~ " [h + 

The energy equation is simply 

l(+) = i(-) (6.20) 

The continuity is simply 

*(+)=* (-) ( 6 - 21 ) 

and the condition of continuity of the stream surface S 2 is given by the <p coordinate of S 2 

r (+) = p(-) ( 6 - 22 ) 

The computer program of reference 64 was modified by working equations derived from the 
preceding basic equations. The modified program was used to carry out the three levels of calculation for 
the DFVLR rotor (ref. 89), for which Laser Two Focus measured internal flow data is available 
(ref. 90). In addition, the experimental Mach number contour plots at several spanwise positions given in 
reference 34 for the peak efficiency operating condition at design speed, the meridional projection of the 
Mach number contours on the S 2 m surface are constructed and shown in figure 6.32(a). The increase of 


+ 1 ! + 1 ln 2l - )2 - 1 (6.19) 

7 + 1 7-1 [7 + I 7 + lJ 
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V_r and s from inlet to outlet, determined from the measured values given in references 34, is used in 
9 

the calculation. 

Figures 6.30 and 6.31 show input values of r, V^r, and s used in the Level II and III calculations. 
From the results obtained in the three levels of calculations, it is immediately seen from figures 6.32 and 
6.33 that the streamlines and Mach number contours obtained in the Level I calculation are quite 
different from Level II and III calculations as well as from the measured result. Mach number contours 
obtained in the Level III calculation agrees best with the measured data. 

It is also seen that the difference between the results obtained in Level II and III calculations are not 
large. The differences in the Sj surface shape and Sj filament thickness, needed for S x calculations, 
are shown in figures 6.33 and 6.34, respectively. 

The relative flow angle obtained in the calculation is compared to that of the mean camber line and the 
measured value (ref. 90) in figure 6.35. The agreement with measured value is quite good. The 
calculated pressure rise is shown in figure 6.36. 

6.7 Direct-Problem Solution of Transonic Flow Along S 2 Stream Filament 

In the case of a direct problem of transonic flow along the S 2 stream filament, it is desirable to put 
the nonorthogonal curvilinear coordinates, x 1 and x 2 on the mid-S 2 surface of the filaments, and the 
x 3 -coordinate perpendicular to the surface (refs. 91 and 92) (fig. 6.37). Then W 3 = 0 If dx 3 = 1, then 
= r n = r B cos (n,p) and 
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The continuity equation, the dynamic equation and stream function velocity, and the principal equations 
are, respectively (ref. 91). 


\[&22 sin tf 12 Wl ) + 
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where « 3 = u cos(n,z) / r n , 
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This principal equation is similar to that on the stream surface, but it contains only two unknowns, 
$ and p. The traditional V-p iteration is used in the solution. 

The method of artificial compressibility is used in a manner quite similar to that in transonic Sj 
flow (section 6.2). In the calculation for the 8:1 presssure ratio centrifugal compressor (fig. 6.38) of 
reference 93, it was found that the viscous effect is large and that only after appropriate entropy 
correction and use of a blockage coefficient of 0.18, about the same value used in reference 93, the 
calculated Mach number at the casing was brought up to a level close to the experimental value (fig. 
6.39). Results of computation for a mass flow of 0.909 kg/s at tip Mach number of 1.25 are given in 
figures 6.40 and 6.41. The increase in V^r along the casing is shown in figure 6.40. The increase in 
across the shock is a small portion of the total increase. Figure 6.41 shows the meridional projection of 
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Mach number contours on S 2 m . It is seen that flow in an inducer of a centrifugal compressor is quite 
similar to that in a rotor of a transonic axial-flow compressor. 
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CHAPTER 7 


THREE-DIMENSIONAL FLOW IN TRANSONIC TURBOMACHINES 


The determination of three-dimensional flow in transonic turbomachines through the use of S x and 
S, stream filaments proceeds in the same manner as in the case of subsonic turbomachines. Various 

techniques described in Chapter 6 for solution of transonic S 1 and S 2 flow are selected for use in the 

three-dimensional solution. Because the blading has a stronger influence than the hub and casing walls 
on the passage shock, the position and shape of the shock are determined in the S x calculation. 

Together with the geometry of the S 2 stream filament and the discontinuous changes of filament 
thickness r, angular momentum V tf r, and entropy s, the position and shape of the shock are carried 
over to the S 2 calculation in the next cycle. This procedure not only was proven to be practical in the 
Sj/S 2 iteration process, but also provides a possibility to modify the geometry of the blade and 
meridional passage for minimizing the loss caused by passage shock. So the three-dimensional blade 
design procedure suggested in this Chapter is even more desirable in the case of transonic turbomachmes. 

7.1 Quasi-Three-Dimensional Flow Field in the DFVLR Rotor 

Quasi-three-dimensional solution refers to the solution obtained by using a number of S t stream 

filaments of revolution and a mid-channel S 2 filament. In reference 95, this approximate method was 

applied to compute the flow field in the DFVLR rotor (ref. 89) operating at the design point, of which 
measured flow field (ref. 90) is available for comparison. The design stagnation pressure ratio is 1.51 and 
the inlet relative Mach number at the blade tip is about 1.37. The V,r distribution at the outlet at 
design point is determined from experimental data and is used as the fixed boundary value of the direct 

problem. 
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For S 2m solution a grid system of 61 stations along x 1 and 11 stations along x 2 is selected on 
the meridional plane. Seven Sj surfaces located at 0, 10, 30, 50, 70, 90, and 100 percent, respectively, of 
the blade height at the inlet station are used in the calculation. The grid system selected on the Sj 
surface has 61 stations along x 1 coordinate and 15 stations along x 2 coordinate. The methods of 
references 42 and 64 are used for Sj and S 2 calculations, respectively. In order to see how high of a 
convergence level can be obtained for this moderately high transonic flow, nine cycles of calculation were 
carried out and the convergence level obtained (relative change of flow variable less than 1 percent) is as 
good as that previously obtained for subsonic flow. The significant advantage of the transonic three- 
dimensional calculation method is that the computation can be carried out at low cost in a reasonable 
time period on a modern microcomputer, which is readily available to all design engineers. For instance, 
the computing time for the transonic S x and S 2 solution is about 30 sec each on a 7-MIP Sun 4 
machine. Thus, the total time required is less than 1 hr. This time decreases in proportion to the 
number of CPU’s in a multiprocessor. If, say, eight Sj flows can be done at the same time on an 8-CPU 
computer, the total time required will be reduced by a factor of eight. Therefore, this relatively 
inexpensive, approximate three-dimensional code is quite suitable for design investigation and trade-off 
studies. The final choice of design configuration may then be checked and refined by a full-three- 
dimensional analysis code or even by a viscous code, which has recently been made available. 

Convergence of Iterative Three-Dimensional Solution 

The key issue for the success of the present method of solving the three-dimensional transonic flow 
by iterative calculation between an S 2 m surface calculation and a number of Sj surface calculations is 
the convergence problem. In the case of a high subsonic flow through a compressor stator (ref. 31) it was 
found that when the entropy values of the two surfaces at the same intersecting point are kept the same, 
a convergent solution with 
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< 2 percent 


was obtained after only four cycles of iteration. For transonic flow it is expected that the convergence of 
solution would be more difficult to obtain than that of the passage shock. However, it is believed that, in 
the present calculation procedure, that the shock discontinuity is determined in the S x solution and that 
the abrupt changes in entropy r and V^r across the shock obtained in the S x solution are taken as 
the input values to the S 2 m solution will help achieve convergence of the three-dimensional iterative 
solution. 


In the example cited above, the criterion |Ar/r| max and |Ar/r| max are used to judge convergence, 
whereas | (M Sj - M s J/M Si J^ is used to judge consistency. Between convergence and consistency, the 
former is more essential. If convergence is not achieved, good consistency is impossible. On the other 
hand, even in the case of a high degree of convergence, there may be an irreducible inconsistency of the 
Mach number. It is believed that consistency of entropy and mass flow coefficient in the Sj and S 2 
calculations is especially important for consistency between the Sj and S 2 solutions. 


Using a proper value of relaxation factor is very important for convergence. During the calculation 
for flow in transonic rotor, the thickness value of the Sj stream filament near the casing, fluctuated and 
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diverged when a relaxation factor was not used (see fig. 7.1). This kind of divergence is something like 
that in the iterative calculation of the following algebraic equation. 


X<»> = f (,<->)) 


when f'(x) <— 1. Hence, in order to obtain and speed up convergence, a small relaxation factor was 
used. 


In the numerical example to be presented in the following section, the same grid pattern in the x 1 
direction is employed for both Sj and S 2 calculations, and the entropy and mass flow coefficients are 
kept consistent in both calculations. It is found that using the relaxation factor of a relatively small 
value is very effective in obtaining quick convergence. The consistency is, at the same time, relatively 
good. In the case of high subsonic flow through a compressor rotor, the relative change in S x filament 
thickness is reduced to below 2 percent and the relative differences between Sj and S 2 Mach numbers 
and flow angles are reduced to 0.9 percent and 0.3 degrees, respectively. In the case of transonic flow in a 
compressor rotor with a tip Mach number of 1.37, the relative change in S x filament thickness is 
reduced to below 1 percent and the relative difference between S x and S 2 Mach numbers at the same 
grid point is below 4 percent. These values of consistency are probably due to the fact that artificial 
density is employed in S x but not in S 2 calculations and that the conservation form of the principal 
equation is employed in Sj but not in S 2 calculations. 
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Variation of Flow Variables During the Convergence Process 


During the first few cycles of calculation there was a general tendency toward convergence. 

However, in the region near the casing, the values of filament thickness soon began to fluctuate. 

Also its magnitude tends to increase (see fig. 7.2). A small relaxation factor, less than 1, was used to 
make the process convergent. It was noticed in the 7th cycle that a small input to Sj calculation 
(the lower solid line in fig. 7.2) yields a large value of output from the S 2 calculation (the upper dashed 
line in fig. 7.2) and vice versa in the 8th cycle. The solution lies between the two. The middle value was 
then taken as input to the 9th cycle and the output obtained is less than 1 percent from the final value. 

The difference in streamlines obtained on the S 2 m surface in the 1st cycle and 9th cycle shown in 
figure 7.3 is seen to lie in the region downstream of the shock. The radius coordinate obtained in the 9th 
cycle in that region is considerably higher than that in the 1st cycle. 

The stream- wise variation of r for the seven Sj surfaces assumed in the 1st cycle and obtained in 
the 9th cycle is shown in figure 7.4, and the value of r and r, respectively, of the streamline at 70- and 
100- percent blade height in the 7th to 9th cycles of calculation is listed, respectively, in tables I and II. 
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TABLE I.— THE VALUES OF THICKNESS (cm) OF THE Sj STREAM FILAMENT AT THE BLADE 

TIP IN THE LAST THREE CYCLES 
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The converging process of V^r is similar to that of r of S r The variation of V^r, relative to its 
value at the exit, along the seven streamlines in the 9th cycle of calculation is compared to that in the 1st 
cycle of calculation in figure 7.5. 

On the other hand, the change in the angular thickness of stream filament S 2 m is relatively simple. 
The stream-wise value on the seven streamlines used as input value to the 9th cycle is shown in 
figure 7.6(b) and is compared to those assumed in the 1st cycle (fig. 7.6(a)). The abrupt change in r 
(and V tf r streamline) is clearly seen in figure 7.5(b) (and 7.3(b), 7.4(b)). 

Figure 7.7(a) shows the Mach number variation along the mean streamlines obtained on the seven 
Sj surfaces in the 9th cycle, whereas figure 7.7(b) shows those along the streamlines on the S 2m surface 
obtained in the 9th cycle. They are very close; the maximum difference between the two being less than 
4 percent. This is the same value previously obtained in the iteration for a subsonic rotor. This means 
that in the present solution of moderately high transonic flow, the convergence and consistency level 
obtained is as good as those in the subsonic case. 

Three-Dimensional Flow Field 

In the approximate three-dimensional solution presented herein, the three-dimensional flow field is 
obtained on the seven S 1 filaments of revolution extending over the whole blade passage (the S 2 m 
calculation is used mainly for obtaining the geometry of the Sj filaments). The constant Mach number 
contours on these S, surfaces can be readily obtained from the calculation. However, in order to make a 
comparison to test data, the contours on S t surfaces at 18, 45, 68, and 89 percent span are obtained 
through interpolation. They are shown in figures 7.8(a) to (d). They clearly show the strong effect of 
the passage shock. The solution of S x surface flow can also be used to obtain the radial variation of 
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flow variables at a number of stations along the flow direction. Two of them, the relative Mach number 
and flow angle j3 are shown in figure 7.9. 

Although only the flow variation on the central S 2 surface is directly available in the approximate 
three-dimensional solution, those on the other S 2 surfaces can be constructed by the data on the 
corresponding streamlines on different Sj surfaces. Constant Mach number contours on the suction and 
pressure surfaces are obtained in this way and are shown together with that on the central S 2 surface in 
figure 7.10. A second passage shock is seen on the suction surface and the central S 2 surface, extending 
only a short distant inward from the casing. 

From data on these S 2 surfaces, stream-wise variations of flow variables on certain streamlines can 
be readily obtained. For instance, the variation of flow angle along three streamlines on the S 2 m 
surface at 18-, 68-, and 89- percent blade height (at inlet) are shown in figure 7.11. 

2 3 

The variation of flow variables can also be shown on spanwise surfaces formed by x and x 
coordinates. The variations of M, and p on these surfaces are particularly useful for observing the 
influence of the passage shock and the magnitude of vorticity or secondary flow. Mach number contours 
on three such surfaces are shown in figure 7.12. On the 10-percent chord surface, the intersection of 
shock surface and this surface is seen to extend from the casing all the way inward to a point very near 
the hub. At the 90-percent chord surface, the flow is almost entirely subsonic. 

In order to more clearly see the effect of the three-dimensional passage shock on the three- 
dimensional flow field, the three-dimensional passage shock is constructed in scale and is shown in 
figure 7.13. At the tip Mach number of 1.37 and hub tip ratio of 0.5 the passage shock extends from 
casing all the way inward almost to the hub. 
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Comparison Between Computed Flow Field with DFVLR Measured Data 


The quasi-three-dimensional transonic solution obtained by the present method is compared with 
available DFVLR measured data in figures 7.8, 7.9 and 7.11. 

The character of the flow field obtained by the present method is seen to be in good agreement with 
that constructed from DFVLR measured data. In the region near 18-percent span, a supersonic zone 
appears near the leading edge of the suction surface. In the middle of the span, a bow wave appears a 
short distance in front of the leading edge and the wave extends to the flow passage to form a passage 
shock. At about two-thirds of the span, the shock appears to be nearly attached to the leading edge. 
Located in the blade tip region, is a second shock downstream of the bow wave/passage shock and the 
shock is nearly normal to the suction surface. The position of the passage shock obtained in the 
calculation is slightly in front of the measured position. It seems that this difference is inherent in the 
solution of transonic flow by the use of the stream function/artificial density method (refs. 41, 42, and 
79). Also the computed shock is more nearly normal to the suction surface, whereas the measured shock 

is more inclined to the suction surface. 


The result calculated by a different method (ref. 34) is included in figure 7.8 for comparison. The 
character of flow is similar, but there is some difference in the magnitude between the two calculated 

results. 

The relatively good agreement between the calculated and measured values of Mach number and 
flow angle in the radial direction (fig. 7.9) and flow angle along the stream-wise direction (fig. 7.11) 
indicate that flow field in the central portion of the three-dimensional passage determined in the present 
method is sufficiently accurate for engineering computation. It seems that a more accurate determination 
of the shock position, such as the separate calculation/shock fitting method (section 6.1), appropriate 
correction for blockage effect at higher radius (section 6.3,) and without assuming axially symmetrical Sj 
stream surface may improve the solutions. 


7.2 Full- Three-Dimensional Transonic Flow in Cas Rotor 


Solution Obtained by Using a Number of Sj and S 2 Stream Filaments 

The method described in Chapter 5 for obtaining three-dimensional flow using a number of S x and 
S 2 stream filaments was employed in reference 96 for calculating full-three-dimensional transonic flow in 
the CAS rotor. A full-three-dimensional calculation was made after a quasi-three-dimensional solution 
was been obtained. Thus, full-three-dimensional calculation can be used as a full-three-dimensional 
analysis solution after the rotor is designed by the quasi-three-dimensional procedure. 

The design parameters of the CAS compressor and details of testing were given in section 5.9. The 
meridional projection of the rotor blades, design streamlines, L2F measuring stations, and the overall 
performance of the rotor were given in Figures 5.42 and 5.43. The radial distribution of stagnation 
pressure and absolute flow angle measured by optical and non-optical devices are shown in figure 5.44. 
Three-dimensional flow calculation has been carried out for 80 and 90 percent of design speed and 
compared to the L2F data available at these two speeds. 

Quasi-Three-Dimensional Solution for Transonic Flow in CAS Rotor 

Quasi-three-dimensional solutions of transonic rotor were obtained for 80- and 90-percent design 
speed. The stagnation pressure and the absolute flow angle, measured by Setra capacitance transistor 
and L2F velocimeter, respectively, are taken as the outlet boundary conditions. The other main input 
data are: flow path in meridional plane, blade geometry (coordinates given on S x surface or on 
manufacturing templates), inlet pressure and temperature, rotor speed, mass flow and mass flow 
coefficient. 
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The First step is to solve the flow along an S 2>m surface about midway between two adjacent 
blades. The computer code of reference 64, with a given absolute outlet flow angle and rotor efficiency, is 
used for this purpose. The isentropic efficiency of the rotor is determined by the torque measurement. 

The radial distribution of rotor efficiency is calculated with the stagnation pressure and temperature 
measurements. With the absolute velocity and flow angle measured by L2F velocimeter, the angular 
momentum V^r is calculated. During calculation, the outlet flow angle a and radial distribution of 
efficiency t) t are updated successively, until V 0 r converges, and the calculated pressure p° is in 
agreement with the experimental values. This process takes only three to five cycles. The set of radial 
distribution values of p°, a, and V,r at the outlet station and the rotor efficiency are used for 
three-dimensional calculation of internal flow field and are compared with the L2F measured internal 
flow field. Figure 7.14 shows the calculated values as compared to the measured values at the outlet 

station. 

From this solution of S 2 m flow, the geometry of the 11 S x surfaces of revolution and the variation of 
the S t stream filaments are obtained. The coordinates of the corresponding blade surface are calculated 
by interpolation. The computer code of reference 41 is used to calculate the S x transonic flow. After 
the solution of the flow on the 11 S x surfaces are obtained, a new S 2 m surface along with the 
corresponding variations of thickness of the S, filament Number 10 in figure 7.15. The angular 

momentum V,r along streamline. Number 2 and Number 10 on the S 2m surface during iteration are 
shown in figure 7.16. In figure 7.17 the relative Mach numbers along the intersecting streamlines 
obtained on g; surface and on S 2 m surface are shown. The closeness between the values shows the 
degree of consistency as well as convergence reached in the solution. 


Table III below gives the geometry and major flow parameters on six S 1 stream filaments obtained 
in the solution 
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TABLE III. — MAJOR PARAMETERS ON SIX S, STREAM FILAMENTS OF REVOLUTION 


s 

M 

degree 

£31 

degree 

* 

Coordinates r,m 

r/r a ' 

r l,0 

degree 

m 

Hi 

m 

m 

m 

HI 

r 0 

11 

1.180 

-66.2 

-66.2 


.3390 

.3354 

1.000 

1.014 

'in- '""M 

.943 

-7.98 

a 

1.089 

-64.2 

-51.2 

.3133 

.3102 

.3078 

1.092 

1.021 


.896 

-5.49 

Bfl 

.993 

-61.3 

-45.7 

.2807 

.2801 

.2788 

1.234 

1.121 

1.071 

.868 

-1.73 

H 

.883 

-57.4 

-38.5 

.2434 

.2460 

.2465 

1.404 

1.305 

1.206 

.859 

2.61 


.758 

-51.2 

-27.6 

.1986 

.2055 

.2096 

1.730 

1.589 

1.433 

.929 

8.62 

H 

.609 

-38.7 

3.8 

.1370 

.1499 

.1628 

2.638 

2.461 

1.929 

.731 

18.41 


V 0 = 1.141 cm. 


In the table values of r 0 , r 0 5 , r^, and r/r 0 show clearly the quasi- three-dimensional assumption of 
the Sj stream surface. Although the S 1 surface is axially-symmetric, the flow on varies in the 
circumferential direction. 


The meridional projection of streamlines on S 2 m is shown in figure 7.18(a). 

The stream-wise variation of the S 2 filament thickness and the angular momentum V^r at the mean 
streamline on six Sj stream surfaces is shown in figures 7.19 and 7.20, respectively. On the Sj surface 
near the hub, the flow is subsonic and these curves resemble those previously obtained in the subsonic 
compressor (for example, fig. 5.40 and 5.39). But for the Sj surface located in the region where a 
passage shock exists, the curves are greatly affected by the shock. 

The circumferential variations of relative Mach number at six axial stations on Five Sj surfaces are 
shown in figure 7.21. 
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Full-Three-Dimensional Flow in CAS Rotor 


A s mentioned in Section 5.7 the full-three-dimensional solution can be made by going into the multi- 
S, calculation immediately after the completion of the multi-Sj calculation in the first cycle interaction 
(A-B-C in fig. 7.22). But, in the present calculation, multi-Sj calculations are carried out after the 
completion of the quasi-three-dimensional solution (A-B-A-B in fig. 7.22). It is found that the iterative 
calculation between six (in the first four cycles) to 11 (in the last two cycles) Sj surfaces and seven S 2 
surfaces reaches an acceptable level of convergence and consistency after six cycles of iteration. There are 
52 (21 in the blade region) and 11 stations, respectively, along x 1 and x 2 coordinates on each Sj 
surface, and 52 and 11 stations, respectively, along x 1 and x 2 coordinates on each S 2 surface. The 
solution of the flow on the seven S 2 surfaces employed in the iterative calculation is based on the 
“inverse type” solution, i.e., V^r obtained in the Sj solution is taken as the input value to the S 2 
solution. The S 2 surface adjacent to the blade surface is placed a short distance from it, and the mass 
flow between them is 5 percent of the total flow. The Sj and S 2 computer codes are connected 
through two interpolation codes so that V 0 i distribution and S 2 filament thickness obtained in the 
solution of flow on Sj surfaces are interpolated for input values to different S 2 surfaces. Similarly the 
geometry of the Sj surfaces and the filament thickness r obtained in the solution of flow on S 2 
surfaces, are interpolated for input values to the S ^ surfaces. 

Similar to the quasi-three-dimensional solution, entropy increases across the rotor at different radii 
are calculated according to the rotor efficiency at these radii. The variation of entropy increased along 
the x 1 direction is obtained by considering that (1) there is no entropy increase in the flow up to shock 
line, (2) entropy increases abruptly across the shock, and (3) from there on it increases linearly to the 
outlet value. In order to be consistent, the value of entropy used in the Sj calculation is the same as 
that in the S 2 calculation at the same grid point. Because of the lack of data on how entropy varies 
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from blade to blade, the circumferential variation of entropy is considered neither on the S ^ surface nor 
on different S 2 surfaces. 

Converging Process 

Since the main difference between a quasi-three-dimensional and a full-three-dimensional solution is 
whether the radius of the Sj surface is or is not circumferentially constant, the convergence in the 
geometry of S x surface is an indication of the convergence of the full-three-dimensional solution. 

Table IV lists the changes in the radii of three Sj surfaces during iteration. The corresponding changes 
in S 1 filament are listed in Table V. 


TABLE IV,— CHANGES IN RADIUS (m) OF GENERAL S 1 SURFACES DURING ITERATION 


Si 

S s 

Station 

cycle 

I 

IV 

VII 


21 

31 

41 

21 

31 

41 

21 

31 

41 

9 

1 

.3111 


mm 

.3125 

.3102 

.3077 


.3098 

.3070 


5 

.3118 


m 

.3131 

.3108 

.3082 


.3108 



6 

.3118 

.3101 

■sal 

.3131 

.3108 

.3082 

.3137 

.3108 

■ 

D 

1 





1 


.2814 

.2802 

.2792 

■I 

5 

7H 




■ 

■ 

.2823 

.2806 

.2793 

■ 

6 

.2791 

9831 

■ 

K3i9 

fl" 

.2794 

.2823 

.2806 



1 



.2473 



.2481 



Hi 

1 

5 

.2430 

.2453 

.2469 

[ .2445 

.2468 

.2476 

11 


mSEM 


6 

.2429 

.2453 

.2469 

.2445 

.2467 

.2476 

11 


■ 
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TABLE V.— CHANGES IN FILAMENT THICKNESS (cm) OF GENERAL 3, SURFACES 


Si 

S a 

I 

IV 

VII 


Station 

cycle 

21 

31 

41 

21 

31 

41 

21 

31 

41 

9 

1 

wmm 

B 

.986 

1.096 

1.035 

.987 

1.087 

1.050 

.987 


5 

UIB 

■31 

.986 


1.039 

.977 

1.061 

1.037 

.955 


6 

B « 

mm 

.985 




1.062 

1.036 

.955 

B 

1 

mgm 

1.131 

1.063 

■ 



mm 

1.059 

.989 

B 

5 

■ 

1.136 

1.075 




BBS 

1.103 

1.034 

B 

6 

1.202 

1.136 

1.075 

■ 

1.121 

1.062 

WSSM 


1.035 


1 

1.388 

1.328 

■W 

B 




■Milk 

mm 


5 

1.366 

1.298 


■ 

■ 


bqeEFSh 


BBS 

L 

6 

1.367 

1.299 

■89 1 






B 


= 1.41 cm. 


Convergence in the radius of the general Sj stream surfaces and in the thickness of the general S ^ 
stream filament is clearly indicated in Tables IV and V. 

The Mach numbers along intersecting lines between S 2m surfaces and Sj surfaces are shown in 
figure 7.17(b). 


In figure 7.23 the Mach numbers along intersecting streamlines between S x surfaces and S 2 surfaces in 
this solution are shown. A comparison between figures 7.17(a) and 7.17(b) and figure 7.23 itself indicates 
that much better convergence and consistency are achieved in the full-three-dimensional solution than 
that achieved in the quasi-three-dimensional solution. 


Geometry of Sj and S 2 Stream Filaments 

Data in tables VI and VII show that the radius of the S 1 surface increases in the circumferential 
direction at the blade leading edge, and the difference decreases toward the blade trailing edge. On the 
other hand, the filament thickness at the leading edge is larger near the pressure surface at large radius 
and larger near the suction surface at small radius. The difference decreases toward the trailing edge. 
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This configuration of the S 1 stream Filaments is a combined result of the flow vorticity and the passage 
shock. In subsonic flow the relative circumferential difference in the radius is found to be about 
0.2 percent, but in transonic flow the maximum value is about 1.3 percent, an order of magnitude higher. 
The geometry of the S 2 stream filaments is given in table VIII. 


TABLE VI.— COORDINATE r (m) OF S l STREAM SURFACES 


Sl, S 2 


Number = 

7 


Number = 4 

station 

I 

II 

III 

IV 

V 

I 

11 

III 

IV 

V 


.2792 

.2802 

.2806 

.2811 

.2815 

.2201 

.2205 

.2214 

.2222 

.2229 


.2792 

.2800 

.2811 

.2820 

.2825 

.2229 

.2230 

.2239 

.2248 

.2259 


.2792 

.2800 

.2806 

.2812 

.2819 

.2243 

.2248 

.2256 

.2263 

.2268 


.2795 

.2802 

.2803 

.2804 

.2802 

.2271 

.2278 

.2284 

.2287 

.2290 

41 

.2789 

.2793 

.2793 

.2793 

.2794 

.2294 

.2299 

.2301 

.2303 

.2307 

45 

.2784 

.2785 

.2785 

.2785 

.2787 

.2304 

.2305 

.2306 

.2305 

.2308 


TABLE VII.— THICKNESS OF S x STREAM FILAMENTS (r/rj 


Sl. s i 

Number 11 (Filament of Revolution) 

Number 7 (General Filament) 

station 

I 

ii 

III 

IV 

V 

I 

II 

III 

IV 

V 


■ESS 

1.017 

1.008 

1.001 

.995 

1.251 

KM 

mi 

ms 

ms 


SIM 

1.043 

1.021 

1.007 

.996 


ms 

ms 

■SOS 

mi 


■£9 

.986 

.970 

1.957 

.957 

1.138 

ms 


mi 

1.104 

41 

mil 

.946 

.954 

.957 

.968 

1.075 

ms 

mi 

1.058 

1.052 

46 

m 

.928 

.929 

.929 

.928 

1.042 

mi 

ms 

1.037 

1.035 



Number 4 (General Filament) 


Number 1 (Filament of Revolution) 


i 

II 

III 

IV 

V 

I 

n 

in 

IV 

IBOll 

16 

ng a 

MEM 


MEM 

1,647 

3.497 

3.488 

3.542 

3.592 


21 

mjm 

mi 


mi 

1.543 

2.876 

2.872 

2.955 

3.053 


31 

warn 

mi 


ms 

1.121 

2.613 

2.603 

2.662 

2.709 

2.738 

41 

mm 

ml 


mEm 

1.297 

1.962 

2.055 

2.099 

2.133 

2.191 

46 

WSm 

1.256 


mi 

1.251 

1.692 

1.703 

1.711 

1.720 

1.741 


r 0 = 1.41 cm. 


Figures 7.24 and 7.25 show, respectively, the 10 S x stream filaments formed by the 11 S 1 stream 
surfaces and the intersecting lines of the 11 Sj stream surfaces with three span-wise surfaces. It is found 
during calculation that the twist of the Sj surface in the radial direction is relatively large at the blade 
leading edge and relatively small at blade trailing edge. For instance, the radii of surface Number 5 
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at the blade leading edge are, respectively. 0.2429 m at the pressure surface and 0.2465 m at the suction 
surface, a difference of 0 to 0.0036 m or 1.48 percent, but at the trailing edge they are, respectively, 

0.2469 and 0.2479, a difference of only 0.001 m or 0.4 percent. 

The geometry of S x stream surfaces is one of the most important three-dimensional flow variables. 
The blade section is now operating inside this filament. The flow characteristics are greatly influenced by 
the variation of filament thickness in stream-wise as well as in circumferential direction. 

Three-Dimensional Flow Field 

The three-dimensional geometry of the S x stream filaments, given in tables VI and VII, and the 
blade sections lying in the filament determine the flow along these filaments and consequently the flow 
over the whole channels (see fig. 5.27). The meridional projection of the streamlines obtained on the 
S„ surface S, filament thickness and V-,r along the mean streamline on six S, stream surfaces are 

2,m & v 

given in figures 7.18(b), 7.19(b) and 7.20(b), respectively. The difference between the solution obtained 
in the present full-three-dimensional solution and the previous quasi-three-dimensional solution is quite 
apparent. The major geometry and inlet and outlet flow conditions of the six S x stream filaments are 
listed in table VI in which a and b denote the two corners at the upstream station, e and f at mid- 
point of chord length and c and d at the downstream station. r c /r a , r^/r^ and f c /r a , T &h y, 
contraction in r and filament thickness r, respectively, on the two sides. 

The flow determined on the six S x filaments is shown in figure 7.26 in terms of constant Mach number 
contours. They are drawn for flow fields including four stations upstream and three stations downstream 
of the blade. The inlet Mach number varies from 0.6 at the hub to 1.25 at the tip. The supersonic 
region begins to appear on surface (3). On surface (7) a bow wave appears and extends to the suction 
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surface as a passage shock. On the tip surface (11) the passage shock hits the suction surface at about 
50- percent chord point. 
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The geometry of the S 2 stream filaments is given in table VIII. Similar to the Sj stream 
filaments, the three-dimensional geometry of the S 2 stream filaments and the meridional wall shape 
determine the flow along these filaments and consequently the flow over the whole channel (see fig. 5.29). 
figure 7.27 shows the constant Mach number contours on the S 2 filaments. Notice that, whereas in the 
quasi-three-dimensional solution only flow on one S 2 surface (the S 2 m surface), which took part in the 
iteration, is directly obtained in the solution, in full-three-dimensional solution the flow on all S 2 
surfaces, which took part in the iteration, are directly obtained in the solution. The circumferentially 
uniform flow at the inlet station becomes nonuniform near the blade leading edge, where tip Mach 
number varies from 1.18 on S 2 (I) to 1.21 on S 2 (VII). It becomes uniform again at three stations 
downstream of the blade. From the pressure surface of one blade toward the suction surface of the 
adjacent blade, the intersection line of the shock and the surface becomes longer and moves downstream. 
The flow condition on S 2 (I) and (VII) should be particularly useful in the investigation of boundary 
layer flow along the blade surfaces. 

U sing data obtained from these two families of stream surfaces a three-dimensional view of the shock 
can easily be constructed (fig. 7.28). The passage shock is seen to extend from the casing part way 
(about two-thirds) down to the hub. The shock is nearly normal to the suction surface and more inclined 
at lower radius. 

Using data obtained on the Sj and S 2 surfaces, constant Mach number contours can also be 
constructed on a number of a third family of surfaces, formed by the and coordinates, nearly 
normal to S| and S 2 surfaces at stations 1 (L.E.), 5, 9, 11, and 21 (T.E) (fig. 7.28). In the stream- 
wise direction the shock intersection line becomes shorter, and the circumferential gradient in Mach 
number becomes smaller. Between station 11 and 21 appears a reverse circumferential Mach number 
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gradient in the lower part of the passage. Consequently there is also a reverse circumferential gradient in 
pressure, which will affect the boundary layer along the hub wall. 

Constant a contours are obtained on five such surfaces and are shown in figure 7.30. At station 
J=0 (one station upstream of the blade) there is a significant difference circumferentially in the upper 
half of the flow channel. At the pressure surface side a is +2°, whereas on the suction surface side, a 
b -2° to -6®. This difference in a means that the air b moving outward on the left and moving 
inward on the right. At station 2 the air on the upper left is now moving inward with -8® to -10®. 

Thus there is a twist of Sj surface between station 0 and station 2. It is nearly uniform 
circumferentially at station 6. At station 10 a sharp change of 2® appears in the flow direction over most 
part of the channel. At station 21 the value a is again nearly uniform in the circumferential direction. 

The stream- wise variations of filament thickness and angular momentum V^r of S 2 stream 
filaments (I), (IV), and (VII) are shown, respectively, in figures 7.31 and 7.32. As for the filament 
thickness r, the variation in the subsonic hub region is similar to that of (P-t)/P, except that there is a 
high peak value at the leading edge on the S 3 (I). In the supersonic upper region the variation of S 2 (I) 
and S 2 (VII) are quite different from (P-t)/P. On S 2 (IV) there is still some resemblance, but it has two 
peaks around the leading edge, indicating the effects of bow wave and passage shock. 

The variation of V^r is somewhat similar. In the hub region it is only slightly modified from the 
usual shape in subsonic flow. But in the upper part there is quite a peak at leading edge near the 
pressure surface and an overshoot at the trailing edge near the suction surface. In the tip region there is 
either a two-step rise or fluctuation in the V^r variation. 
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Difference Between Three-Dimensional and Quasi-Three-Dimensional Solution 


The difference between three-dimensional and quasi-three-dimensional solution is entirely due to the 
difference in the geometry of the Sj stream filament. Tables III and VIII show the major difference 
between the Sj stream filaments involved in the two calculations. 

First, there is some difference in the geometry of S 2 m in the radial coordinate and in the 

contraction of the filaments thickness. Figure 7.19 shows the difference in the streamwise variation of the 

S* _ filament thickness. 
z,m 

A rather complicated variation in the full- three-dimensional flow is much simplified in the quasi- 
three-dimensional flow. The streamwise variations of the important flow variable V^r on the S 2 >m 
surface are compared in figure 7.20. A higher overshoot exists in the quasi-three-dimensional solution. 

The streamwise variations of Mach number along the intersecting streamline of the Sj and S 2>ra 
surfaces, obtained in the solution on the two surfaces, are compared in figure 7.17. The full-three- 
dimensional solution has reached a higher level of consistency and convergence than the quasi-three- 
dimensional solution. The former, in general, gives a slightly higher peak Mach number than the latter. 
This difference is also noticeable in the comparison of the constant Mach number contour maps on the 
S 2 m surface shown in figure 7.34. 

Circumferential variations of Mach number at six stations on five Sj surfaces obtained in full-three- 
dimensional and quasi- three-dimensional solution are compared in figure 7.21. In addition to the factors 
influencing the solution mentioned above, the circumferential variation in the S x geometry in the full- 
three-dimensional solution also influences the result. It can be seen from tables V, VI, and VIII that for 
Sj surfaces above Sj(7), the filament thickness and its contraction on the side near the pressure surface 
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is significantly different from those on the side near the suction surface. The combined effect on the 
velocity distribution is that, on Sj surfaces above Sj(7), the Mach number of the full-three-dimensional 
solution is higher in the region from mid-channel to suction surface and lower in the region from mid- 
channel to pressure surface, than that of the quasi-three-dimensional solution. 

Comparison with L2F Measurement 

Efforts were made at IETP to obtain internal flow on Sj surface (4), and (7) to (10) by the use of 
a Laser-2-Focus (L2F) velocimeter developed by DFVLR. Measured values are shown in figure 7.21. In 
general the measured values are a little higher than the calculated values. Constant Mach number 
contours, constructed from the L2F measured data, are shown in figure 7.33. The trend is similar, but 
again the measured Mach number is a little higher than the calculated value. 

7.3 Some Remarks 


The general theory of three-dimensional flow in subsonic and supersonic turbomachinery based on 
the iterative solution between Sj and S 2 stream filaments has been successfully extended to a transonic 
flow regime. In this extension, some of the recently developed simple and quick methods for solving 
transonic flow along Sj and S 2 stream filaments are utilized. For the CAS transonic rotor, three- 
dimensional solution is obtained after six cycles of iteration by using the quasi-three-dimensional solution 
as the starting value. It has a higher level of convergence and consistency than the quasi-three- 
dimensional solution. The difference between three-dimensional and quasi-three-dimensional solution is: 
The former gives a higher peak Mach number than the latter. This difference is mainly due to the 
neglect of the circumferential variation of the Sj filament thickness in the quasi-three-dimensional 
solution. A comparison between the theoretical solution and the L2F measurement shows that the 
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character of the transonic flow including the three-dimensional shock structure is in fair agreement, and 
the calculated velocity is a little lower than that measured by L2F over most of the flow region, except 
at the 20-percent chord station in the upper half region near the suction surface, where the measured 
velocity is a little lower. 

A practical three-dimensional computer aided design and analysis (three-dimensional-CADA) system 
has been constructed by the use of the general Sj code and the general S 2 code. The recommended 
procedure of three-dimensional-CADA calculation follows. 

(1) Prescribed the through-flow design requirements on the S 2 ni surface and obtain the S 2 m 
solution 

(2) Design the blade section on a number of general Sj stream filaments 

(3) Stack the blade sections 

(4) Obtain the three-dimensional analysis for the blade just designed 

(5) Modify step 1 or step 2 and repeat steps 3 and 4 until the results obtained meet the design 
requirements 

It is believed that the configuration of the three-dimensional geometry of the Sj stream filaments, 
obtained during the iterative calculation between the and S 2 flow, to a large extent, determines the 
performance of the blade section lying in the stream filament. With this information available better 
design of the blade section and, consequently, of the whole blade can be realized. Another advantage of 
the presented method is that the computer storage required is very small and the amount of CPU time is 
also small so it is feasible to carry out all of the computations on a microcomputer. 
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Figure 3.5. — Element of Sg stream filament. 


Figure 3.6— S 2 stream surface and S 2 stream filament 



Figure 3.7.— Orthogonal curvilinear coordinates (, <f, on S, surface of revolution. 
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Figure 4.1 . — Non-orthogonal curvilinear coordinates employed 
In the Investigation of S 1 flow. 
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Figure 4.2. — General coordinates x 1 and x 2 , base vectors e<| 
and e 2 , reciprocal vector* e 1 and e 2 and contravariant and 
covariant component* of vector W. 



Figure 4.4.— Unit vector and vorticity normal to 
Si surface. 






X2 



Figure 4.7— Line integration around grid point E. 


184 



1 




Figure 5.1 .—Mean stream surface and mean 
streamline. 
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Figure 5.2.— Variations of velocity on blade surface and mean Figure 5.3.— Streamline distribution in turbine cascade, 

streamline. 
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Figure 5.6.— Mean streamline in the flow passage. 
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Figure 5.7.— Streamline in a radial-flow impeller. 


Figure 5.8 — Variation of meridional velocity on mean streamline. 
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Figure 5.9. — Comparison of velocity distribution around a 
transonic turbine blade. 
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Figure 5.1 1 . — Family of S 2 surface formed progres- 
sively from S 2 . 



Figure 5.12.— Projection of S 2 m streamlines on meridional plane 
(CAS compressor). 




Figure 5.13.— Blade shape and velocity distribution on J = 3 S 1 
surface. 


Figure 5.1 5.— Transformation of coordinates. 


Figure 5.16.— Projection of S 2 m streamlines on mertdonal 
plane (CAS compressor). 






Figure 5.1 7.— Blade shape and velocity distribution on k-7 
coordinate surface. 
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Figure 5.1 8. — Comparison of mach number distribution on k-9 
S 1 surface. 




Figure 5.20.— Body-fitting quasi-orthogonai C-grid around the 
blade leading edge. 
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Figure 5.31 .—Calculation steps in quasi-3-D blade design. 
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Figure 5.32.— Determination of flow condition far upstream and downstream of the 
embedded blade row. 
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gure 5.33. — Projection of stator blade and streamline 
distribution on meridional plane. 
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Figure 5.34.— Effect of entropy on consistency 
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Figure 5.40. — Variation of B/B k and P/( p-t) on mean streamline 
(mid-span surface). 



Figure 5.41 .—Mach number distribution over blade surface 
(mid-span section). 



Figure 5.42.— Meridional projection of rotor blade, streamlines 
and L2F measuring locations. 
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Figure 5.43.— Compressor overall performance map. 
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Figure 5.46.— Rate of convergence. 
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Figure 5.47.— Radial variation of p. Figure 5.48.— Radial variations of relative Mach number. 
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Figure 5.49. — Radial variation of absolute Mach number at Stations 23 - 28. 
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Figure 5.50,-Axial variation of relative flow angle. of cons,ant relative Mach number on Na 3 
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Figure 6.1 — Bow wave/passage shock and characteristic net in 
a compressor cascade with supersonic inlet velocity. 



Figure 6.2 — Coordinate system used downstream 
of passage shock. 
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(a) Velocity W m/s. (b) Flow direction 0. 

Figure 6.5.— Comparison between result calculated by R-H relations and result obtained in downstream subsonic solution 
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Figure 6.6. — Mach number distribution for turbine cascade. 
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Figure 6.7.— Mach number distribution for T n cascade. 
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Figure 6.8.— Modification Mach number distribution for T 1 
cascade. 
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Figure 6.10. — Pressure distribution around the blade. 
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Figure 6.1 1 .—Flow variation along three streamlines 
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Figure 6.1 2. — Pressure distribution around the blade. Figure 6.1 3. — Flow variation along three streamlines. 
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Figure 6.14. — Pressure distribution around blade. 


205 



Figure 6.1 6. — Variation of r and ^ for a general S 1 surface. 
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Figure 6.21 .—Meridional projection of Mach number Figure 6.22.— Variation of radial thickness of S, stream filament, 

contours on S 2 rTV 



Fjgure 6.23. Variation of Relative flow angle. Figure 6.24.— Variation In relative Mach number. 
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Figure 6.25.— Variation of absolute flow angle. 
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Figure 6.26.— Variation in absolute Mach number. 
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Figure 6.27.— Variation of V 0 r. 
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Figure 6.28.— Variation of pressure ratio. 
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(a) Grid used in calculations I and II, 61x11 
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(b) Grid used In calculation III, 49x11 . 

Figure 6.29. — Grid system and project of blade and shock on 
meridional plane. 
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Figure 6.35. — Comparison of relative flow angle. 


Figure 6.36. — Pressure rise on S 2 t r 
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Figure 6.37.— N on-orthogonal curvilinear coordinate system 
used tor S 2 stream surface. 
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Figure 6.38 — Meridonal projection of mesh on S 
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Figure 7.1 .—Local fluctuation of ^ stream filament thick- 
ness during iteration and the effect of relaxation factor 
(sketch). 
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Figure 7.3.— Meridional projection of sti 
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Figure 7.8.— Mach number contours on S, 


surface: (1 ) DFVLR measured data; (2) calculated by present method; (3) from Ref. 34. 
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Figure 7.1 0— Constant Mach number contours on three S 2 surfaces. 
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Figure 7.24.— Ten S n stream filaments. 
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Figure 7 . 25 .— Intersection of S 1 stream surfaces with three 
span-wise surfaces. 
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(a) Full 3-D. (b) Quasl-3-D. ( c ) L2F - 

Figure 7.33.— Comparison of Mach number contours on S 2m . 
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